RINGMesh  Version 5.0.0
A programming library for geological model meshes
io_tsolid.cpp
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2012-2017, Association Scientifique pour la Geologie et ses
3  * Applications (ASGA). All rights reserved.
4  *
5  * Redistribution and use in source and binary forms, with or without
6  * modification, are permitted provided that the following conditions are met:
7  * * Redistributions of source code must retain the above copyright
8  * notice, this list of conditions and the following disclaimer.
9  * * Redistributions in binary form must reproduce the above copyright
10  * notice, this list of conditions and the following disclaimer in the
11  * documentation and/or other materials provided with the distribution.
12  * * Neither the name of ASGA nor the
13  * names of its contributors may be used to endorse or promote products
14  * derived from this software without specific prior written permission.
15  *
16  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
17  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO,
18  * THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
19  * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL ASGA BE LIABLE FOR ANY DIRECT,
20  * INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
21  * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
22  * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
23  * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
24  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
25  * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26  *
27  * http://www.ring-team.org
28  *
29  * RING Project
30  * Ecole Nationale Superieure de Geologie - GeoRessources
31  * 2 Rue du Doyen Marcel Roubault - TSA 70605
32  * 54518 VANDOEUVRE-LES-NANCY
33  * FRANCE
34  */
35 
36 namespace {
37  using namespace RINGMesh;
38 
39  class TSolidIOHandler final : public GeoModelIOHandler< 3 > {
40  public:
41  TSolidIOHandler()
42  :
43  GeoModelIOHandler< 3 >(),
44  out_(),
45  numeric_like_vertex_attribute_names_(),
46  vertex_attribute_dimensions_(),
47  numeric_like_cell_attribute_names_(),
48  cell_attribute_dimensions_(),
49  vertex_exported_(),
50  atom_exported_(),
51  vertex_exported_id_(),
52  atom_exported_id_()
53  {
54  }
55  void load( const std::string& filename, GeoModel3D& geomodel ) final
56  {
57  std::ifstream input( filename.c_str() );
58  if( !input ) {
59  throw RINGMeshException( "I/O", "Failed loading geomodel from file ",
60  filename );
61  }
62  GeoModelBuilderTSolid builder( geomodel, filename );
63  builder.build_geomodel();
64  }
65  void save( const GeoModel3D& geomodel, const std::string& filename ) final
66  {
67  ringmesh_assert( !out_.is_open() );
68  out_.open( filename.c_str() );
69  ringmesh_assert( out_.is_open() );
70  out_.precision( 16 );
71 
72  fill_top_header( geomodel );
73 
74  const GeoModelMesh3D& mesh = geomodel.mesh;
75  //mesh.set_duplicate_mode( GeoModelMeshCells::ALL ) ;
76 
77  fill_vertex_attribute_header( geomodel );
78  fill_cell_attribute_header( geomodel );
79 
80  vertex_exported_.resize( mesh.vertices.nb(), false );
81  atom_exported_.resize( mesh.cells.nb_duplicated_vertices(), false );
82  vertex_exported_id_.resize( mesh.vertices.nb(), NO_ID );
83  atom_exported_id_.resize( mesh.cells.nb_duplicated_vertices(), NO_ID );
84  for( index_t r = 0; r < geomodel.nb_regions(); r++ ) {
85  const auto& region = geomodel.region( r );
86  export_one_region( region );
87  }
88 
89  export_model( geomodel );
90  export_model_region( geomodel );
91 
92  out_ << "END" << std::endl;
93  }
94  private:
95  const GeoModelMesh3D& geomodel_mesh( const Region3D& region ) const
96  {
97  return region.geomodel().mesh;
98  }
99  void fill_top_header( const GeoModel3D& geomodel )
100  {
101  ringmesh_assert( out_.is_open() );
102  // Print Model3d headers
103  out_ << "GOCAD TSolid 1" << EOL << "HEADER {" << EOL << "name:"
104  << geomodel.name() << EOL << "}" << EOL;
105 
106  out_ << "GOCAD_ORIGINAL_COORDINATE_SYSTEM" << EOL << "NAME Default"
107  << EOL << "AXIS_NAME \"X\" \"Y\" \"Z\"" << EOL
108  << "AXIS_UNIT \"m\" \"m\" \"m\"" << EOL << "ZPOSITIVE Elevation"
109  << EOL << "END_ORIGINAL_COORDINATE_SYSTEM" << EOL;
110  }
111  void fill_vertex_attribute_header( const GeoModel3D& geomodel )
112  {
113  ringmesh_assert( out_.is_open() );
114  std::vector< bool > is_integer_like_attribute;
115  for( index_t reg_i = 0; reg_i < geomodel.nb_regions(); ++reg_i ) {
116  const Region3D& cur_reg = geomodel.region( reg_i );
117  GEO::AttributesManager& reg_vertex_attr_mgr =
118  cur_reg.vertex_attribute_manager();
119  GEO::vector< std::string > att_v_names;
120  reg_vertex_attr_mgr.list_attribute_names( att_v_names );
121  ringmesh_assert( att_v_names.size() == reg_vertex_attr_mgr.nb() );
122  for( const std::string& cur_att_v_name : att_v_names ) {
123 
124  if( cur_att_v_name == "point" ) {
125  continue;
126  }
127 
128  if( contains( numeric_like_vertex_attribute_names_,
129  cur_att_v_name ) ) {
130  continue;
131  }
132 
133  const GEO::AttributeStore* attr_store =
134  reg_vertex_attr_mgr.find_attribute_store( cur_att_v_name );
135  ringmesh_assert( attr_store != nullptr );
136 
137  if( !GEO::ReadOnlyScalarAttributeAdapter::can_be_bound_to(
138  attr_store ) ) {
139  continue;
140  }
141 
142  numeric_like_vertex_attribute_names_.push_back( cur_att_v_name );
143  index_t cur_dim = attr_store->dimension();
144  vertex_attribute_dimensions_.push_back( cur_dim );
145 
146  const GEO::ReadOnlyScalarAttributeAdapter adapter(
147  reg_vertex_attr_mgr, cur_att_v_name );
148  ringmesh_assert(
149  adapter.element_type()
150  != GEO::ReadOnlyScalarAttributeAdapter::ET_NONE );
151  is_integer_like_attribute.push_back(
152  adapter.element_type()
153  < GEO::ReadOnlyScalarAttributeAdapter::ET_FLOAT32 );
154  }
155  }
156 
157  const auto nb_numeric_like_vertex_attribute_names =
158  static_cast< index_t >( numeric_like_vertex_attribute_names_.size() );
159  if( !numeric_like_vertex_attribute_names_.empty() ) {
160  ringmesh_assert( nb_numeric_like_vertex_attribute_names > 0 );
161  out_ << "PROPERTIES";
162  for( const std::string& cur_num_like_v_att_name : numeric_like_vertex_attribute_names_ ) {
163  out_ << " " << cur_num_like_v_att_name;
164  }
165  out_ << EOL;
166  out_ << "PROP_LEGAL_RANGES";
167  for( auto i : range( nb_numeric_like_vertex_attribute_names ) ) {
168  ringmesh_unused( i );
169  out_ << " **none** **none**";
170  }
171  out_ << EOL;
172  out_ << "NO_DATA_VALUES";
173  write_no_data_value( nb_numeric_like_vertex_attribute_names );
174  out_ << EOL;
175  out_ << "READ_ONLY";
176  for( auto i : range( nb_numeric_like_vertex_attribute_names ) ) {
177  ringmesh_unused( i );
178  out_ << " 1";
179  }
180  out_ << EOL;
181  out_ << "PROPERTY_CLASSES";
182  for( const std::string& cur_num_like_v_att_name : numeric_like_vertex_attribute_names_ ) {
183  out_ << " " << cur_num_like_v_att_name;
184  }
185  out_ << EOL;
186  out_ << "PROPERTY_KINDS";
187  for( auto i : range( nb_numeric_like_vertex_attribute_names ) ) {
188  if( is_integer_like_attribute[i] ) {
189  out_ << " \"Number\"";
190  } else {
191  out_ << " \"Real Number\"";
192  }
193  }
194  out_ << EOL;
195  out_ << "PROPERTY_SUBCLASSES";
196  for( auto i : range( nb_numeric_like_vertex_attribute_names ) ) {
197  ringmesh_unused( i );
198  out_ << " QUANTITY Float";
199  }
200  out_ << EOL;
201  out_ << "ESIZES";
202  for( const auto& cur_v_att_dim : vertex_attribute_dimensions_ ) {
203  out_ << " " << GEO::String::to_string( cur_v_att_dim );
204  }
205  out_ << EOL;
206  out_ << "UNITS";
207  for( auto i : range( nb_numeric_like_vertex_attribute_names ) ) {
208  ringmesh_unused( i );
209  out_ << " unitless";
210  }
211  out_ << EOL;
212  for( auto i : range( nb_numeric_like_vertex_attribute_names ) ) {
213  out_ << "PROPERTY_CLASS_HEADER "
214  << numeric_like_vertex_attribute_names_[i] << " {" << EOL;
215  if( is_integer_like_attribute[i] ) {
216  out_ << "kind: Number" << EOL;
217  } else {
218  out_ << "kind: Real Number" << EOL;
219  }
220  out_ << "unit: unitless" << EOL;
221  out_ << "}" << EOL;
222  }
223  }
224  }
225 
226  void fill_cell_attribute_header( const GeoModel3D& geomodel )
227  {
228  ringmesh_assert( out_.is_open() );
229  std::vector< bool > is_integer_like_attribute;
230  for( const auto& cur_reg : region_range<3>( geomodel ) ) {
231  auto& reg_cell_attr_mgr =
232  cur_reg.cell_attribute_manager();
233  GEO::vector< std::string > att_c_names;
234  reg_cell_attr_mgr.list_attribute_names( att_c_names );
235  ringmesh_assert( att_c_names.size() == reg_cell_attr_mgr.nb() );
236  for( const auto& cur_att_c_name : att_c_names ) {
237 
238  if( contains( numeric_like_cell_attribute_names_,
239  cur_att_c_name ) ) {
240  continue;
241  }
242 
243  const auto* attr_store =
244  reg_cell_attr_mgr.find_attribute_store( cur_att_c_name );
245  ringmesh_assert( attr_store != nullptr );
246 
247  if( !GEO::ReadOnlyScalarAttributeAdapter::can_be_bound_to(
248  attr_store ) ) {
249  continue;
250  }
251 
252  numeric_like_cell_attribute_names_.push_back( cur_att_c_name );
253  auto cur_dim = attr_store->dimension();
254  cell_attribute_dimensions_.push_back( cur_dim );
255 
256  const GEO::ReadOnlyScalarAttributeAdapter adapter(
257  reg_cell_attr_mgr, cur_att_c_name );
258  ringmesh_assert(
259  adapter.element_type()
260  != GEO::ReadOnlyScalarAttributeAdapter::ET_NONE );
261  is_integer_like_attribute.push_back(
262  adapter.element_type()
263  < GEO::ReadOnlyScalarAttributeAdapter::ET_FLOAT32 );
264  }
265  }
266 
267  const auto nb_numeric_like_cell_attribute_names =
268  static_cast< index_t > ( numeric_like_cell_attribute_names_.size() );
269  if( !numeric_like_cell_attribute_names_.empty() ) {
270  out_ << "TETRA_PROPERTIES";
271  for( const auto& cur_num_like_c_att_name : numeric_like_cell_attribute_names_ ) {
272  out_ << " " << cur_num_like_c_att_name;
273  }
274  out_ << EOL;
275  out_ << "TETRA_PROP_LEGAL_RANGES";
276  for( auto i : range( nb_numeric_like_cell_attribute_names ) ) {
277  ringmesh_unused( i );
278  out_ << " **none** **none**";
279  }
280  out_ << EOL;
281  out_ << "TETRA_NO_DATA_VALUES";
282  write_no_data_value( nb_numeric_like_cell_attribute_names );
283  out_ << EOL;
284  out_ << "READ_ONLY";
285  for( auto i : range( nb_numeric_like_cell_attribute_names ) ) {
286  ringmesh_unused( i );
287  out_ << " 1";
288  }
289  out_ << EOL;
290  out_ << "TETRA_PROPERTY_CLASSES";
291  for( const auto& cur_num_like_c_att_name : numeric_like_cell_attribute_names_ ) {
292  out_ << " " << cur_num_like_c_att_name;
293  }
294  out_ << EOL;
295  out_ << "TETRA_PROPERTY_KINDS";
296  for( auto i : range( nb_numeric_like_cell_attribute_names ) ) {
297  if( is_integer_like_attribute[i] ) {
298  out_ << " \"Number\"";
299  } else {
300  out_ << " \"Real Number\"";
301  }
302  }
303  out_ << EOL;
304  out_ << "TETRA_PROPERTY_SUBCLASSES";
305  for( auto i : range( nb_numeric_like_cell_attribute_names ) ) {
306  ringmesh_unused( i );
307  out_ << " QUANTITY Float";
308  }
309  out_ << EOL;
310  out_ << "TETRA_ESIZES";
311  for( const auto& cur_cell_attr_dim : cell_attribute_dimensions_ ) {
312  out_ << " " << std::to_string( cur_cell_attr_dim );
313  }
314  out_ << EOL;
315  out_ << "TETRA_UNITS";
316  for( auto i : range( nb_numeric_like_cell_attribute_names ) ) {
317  ringmesh_unused( i );
318  out_ << " unitless";
319  }
320  out_ << EOL;
321  for( auto i : range( nb_numeric_like_cell_attribute_names ) ) {
322  out_ << "TETRA_PROPERTY_CLASS_HEADER "
323  << numeric_like_cell_attribute_names_[i] << " {" << EOL;
324  if( is_integer_like_attribute[i] ) {
325  out_ << "kind: Number" << EOL;
326  } else {
327  out_ << "kind: Real Number" << EOL;
328  }
329  out_ << "unit: unitless" << EOL;
330  out_ << "}" << EOL;
331  }
332  }
333  }
334 
335  void export_one_region( const Region3D& region )
336  {
337  ringmesh_assert( out_.is_open() );
338  out_ << "TVOLUME " << region.name() << EOL;
339  export_region_vertices( region );
340  export_tetrahedra( region );
341  }
342 
343  void export_region_vertices( const Region3D& region )
344  {
345  // Export not duplicated vertices
346  for( auto c : range( region.nb_mesh_elements() ) ) {
347  export_region_cell_vertices( region, c );
348  }
349  }
350 
351  void export_region_cell_vertices( const Region3D& region, index_t c )
352  {
353  const auto& mesh = geomodel_mesh( region );
354  auto cell = mesh.cells.cell( region.gmme().index(), c );
355  auto cell_center = mesh.cells.barycenter( cell );
356  for( auto v : range( mesh.cells.nb_vertices( cell ) ) ) {
357  export_region_cell_vertex( region, cell, v, cell_center );
358  }
359  }
360 
361  void export_region_cell_vertex(
362  const Region3D& region,
363  index_t cell,
364  index_t v,
365  const vec3& cell_center )
366  {
367  ringmesh_assert( out_.is_open() );
368  const auto& mesh = geomodel_mesh( region );
369  auto atom_id = mesh.cells.duplicated_corner_index( { cell, v } );
370  if( atom_id == NO_ID ) {
371  auto vertex_id = mesh.cells.vertex( { cell, v } );
372  if( vertex_exported_[vertex_id] ) return;
373  vertex_exported_[vertex_id] = true;
374  vertex_exported_id_[vertex_id] = nb_vertices_exported_;
375  // PVRTX keyword must be used instead of VRTX keyword because
376  // properties are not read by Gocad if it is VRTX keyword.
377  out_ << "PVRTX " << nb_vertices_exported_++ << " "
378  << mesh.vertices.vertex( vertex_id );
379 
381  export_region_cell_vertex_attributes( region, vertex_id,
382  cell_center );
383  out_ << EOL;
384  }
385  }
386 
387  void export_region_cell_vertex_attributes(
388  const Region3D& region,
389  index_t vertex_id,
390  const vec3& cell_center )
391  {
392  ringmesh_assert( out_.is_open() );
393  auto& reg_vertex_attr_mgr =
394  region.vertex_attribute_manager();
395  auto vertex_id_in_reg = find_gmm_cell_in_gm_region( region, vertex_id,
396  cell_center );
397  ringmesh_assert( vertex_id_in_reg != NO_ID );
398 
399  for( auto attr_dbl_itr : range(
400  static_cast< index_t >( numeric_like_vertex_attribute_names_.size() ) ) ) {
401  if( reg_vertex_attr_mgr.is_defined(
402  numeric_like_vertex_attribute_names_[attr_dbl_itr] ) ) {
403  ringmesh_assert(
404  reg_vertex_attr_mgr.find_attribute_store(
405  numeric_like_vertex_attribute_names_[attr_dbl_itr] )
406  != nullptr );
407  ringmesh_assert(
408  reg_vertex_attr_mgr.find_attribute_store(
409  numeric_like_vertex_attribute_names_[attr_dbl_itr] )->dimension()
410  == vertex_attribute_dimensions_[attr_dbl_itr] );
411  ringmesh_assert(
412  GEO::ReadOnlyScalarAttributeAdapter::can_be_bound_to(
413  reg_vertex_attr_mgr.find_attribute_store(
414  numeric_like_vertex_attribute_names_[attr_dbl_itr] ) ) );
415  GEO::ReadOnlyScalarAttributeAdapter cur_attr(
416  reg_vertex_attr_mgr,
417  numeric_like_vertex_attribute_names_[attr_dbl_itr] );
418  for( auto dim_itr :range(vertex_attribute_dimensions_[attr_dbl_itr] ) ) {
419  out_ << " "
420  << cur_attr[vertex_id_in_reg
421  * vertex_attribute_dimensions_[attr_dbl_itr]
422  + dim_itr];
423  }
424  } else {
425  write_no_data_value(
426  vertex_attribute_dimensions_[attr_dbl_itr] );
427  }
428  }
429  }
430 
431  index_t find_gmm_cell_in_gm_region(
432  const Region3D& region,
433  index_t vertex_id,
434  const vec3& cell_center ) const
435  {
436  const auto& mesh = geomodel_mesh( region );
442  const auto& gme_vertices =
443  mesh.vertices.gme_vertices( vertex_id );
444  for( const auto& cur_gme_vertex : gme_vertices ) {
445  if( cur_gme_vertex.gmme != region.gmme() ) {
446  continue;
447  }
448  auto cells_around_vertex =
449  region.cells_around_vertex( cur_gme_vertex.v_index, NO_ID );
453  for( const auto& cur_cell_around_vertex : cells_around_vertex ) {
454  auto center = region.mesh_element_barycenter(
455  cur_cell_around_vertex );
456  if( ( center - cell_center ).length()
457  < region.geomodel().epsilon() ) {
458  return cur_gme_vertex.v_index;
459  }
460  }
461  }
462  ringmesh_assert_not_reached;
463  return NO_ID;
464  }
465 
466  void export_tetrahedra( const Region3D& region )
467  {
468  // Export duplicated vertices
469  export_duplicated_vertices();
470 
471  // Mark if a boundary is ending in the region
472  mark_boundary_ending_in_region( region );
473 
474  const auto& mesh = geomodel_mesh( region );
475  for( auto c : range(region.nb_mesh_elements() ) ) {
476  out_ << "TETRA";
477  auto cell = mesh.cells.cell( region.gmme().index(), c );
478  export_tetra( region, cell );
479  out_ << EOL;
480  out_ << "# CTETRA " << region.name();
481  export_ctetra( region, c );
482  out_ << EOL;
483  }
484  }
485 
486  void export_duplicated_vertices() const
487  {
488  /*for( index_t c = 0; c < region.nb_mesh_elements(); c++ ) {
489  index_t cell = mesh.cells.cell( r, c ) ;
490  for( index_t v = 0; v < mesh.cells.nb_vertices( cell ); v++ ) {
491  index_t atom_id ;
492  if( mesh.cells.is_corner_duplicated( cell, v, atom_id ) ) {
493  if( atom_exported[atom_id] ) continue ;
494  atom_exported[atom_id] = true ;
495  atom_exported_id[atom_id] = nb_vertices_exported ;
496  index_t vertex_id = mesh.cells.vertex( cell, v ) ;
497  out << "ATOM " << nb_vertices_exported++ << " "
498  << vertex_exported_id[vertex_id] << EOL ;
499  }
500  }
501  }*/
502  }
503 
506  void mark_boundary_ending_in_region( const Region3D& region )
507  {
508  std::map< index_t, index_t > sides;
509  for( auto s : range( region.nb_boundaries() ) ) {
510  if( sides.count( region.boundary_gmme( s ).index() ) > 0 ) {
511  // a surface is encountered twice, it is ending in the region
512  sides[region.boundary_gmme( s ).index()] = 2;
513  } else {
514  sides[region.boundary_gmme( s ).index()] = region.side( s );
515  }
516  }
517  }
518 
519  void export_tetra( const Region3D& region, index_t cell )
520  {
521  export_tetra_coordinates( region, cell );
523  export_tetra_attributes( region, cell );
524  }
525 
526  void export_tetra_coordinates( const Region3D& region, index_t cell )
527  {
528  ringmesh_assert( out_.is_open() );
529  const auto& mesh = geomodel_mesh( region );
530  for( auto v : range(
531  region.geomodel().mesh.cells.nb_vertices( cell ) ) ) {
532  auto atom_id = mesh.cells.duplicated_corner_index(
533  { cell, v } );
534  if( atom_id == NO_ID ) {
535  auto vertex_id = mesh.cells.vertex(
536  { cell, v } );
537  out_ << " " << vertex_exported_id_[vertex_id];
538  } else {
539  out_ << " " << atom_exported_id_[atom_id];
540  }
541  }
542  }
543 
544  void export_tetra_attributes( const Region3D& region, index_t cell )
545  {
546  ringmesh_assert( out_.is_open() );
547  auto& reg_cell_attr_mgr =
548  region.cell_attribute_manager();
549  const auto& mesh = geomodel_mesh( region );
550  auto center = mesh.cells.barycenter( cell );
551  const auto c_in_reg =
552  region.cell_nn_search().get_neighbors( center,
553  region.geomodel().epsilon() );
554  ringmesh_assert( c_in_reg.size() == 1 );
555  for( auto attr_dbl_itr : range(
556  numeric_like_cell_attribute_names_.size() ) ) {
557  if( reg_cell_attr_mgr.is_defined(
558  numeric_like_cell_attribute_names_[attr_dbl_itr] ) ) {
559  ringmesh_assert(
560  reg_cell_attr_mgr.find_attribute_store(
561  numeric_like_cell_attribute_names_[attr_dbl_itr] )
562  != nullptr );
563  ringmesh_assert(
564  reg_cell_attr_mgr.find_attribute_store(
565  numeric_like_cell_attribute_names_[attr_dbl_itr] )->dimension()
566  == cell_attribute_dimensions_[attr_dbl_itr] );
567  ringmesh_assert(
568  GEO::ReadOnlyScalarAttributeAdapter::can_be_bound_to(
569  reg_cell_attr_mgr.find_attribute_store(
570  numeric_like_cell_attribute_names_[attr_dbl_itr] ) ) );
571  GEO::ReadOnlyScalarAttributeAdapter cur_attr( reg_cell_attr_mgr,
572  numeric_like_cell_attribute_names_[attr_dbl_itr] );
573  for( auto dim_itr : range(
574  cell_attribute_dimensions_[attr_dbl_itr] ) ) {
575  out_ << " "
576  << cur_attr[c_in_reg[0]
577  * cell_attribute_dimensions_[attr_dbl_itr] + dim_itr];
578  }
579  } else {
580  write_no_data_value( cell_attribute_dimensions_[attr_dbl_itr] );
581  }
582  }
583  }
584 
585  void export_ctetra( const Region3D& region, index_t c )
586  {
587  ringmesh_assert( out_.is_open() );
588  const auto& mesh = geomodel_mesh( region );
589  for( auto f : range( mesh.cells.nb_facets( c ) ) ) {
590  out_ << " ";
591  index_t polygon { NO_ID };
592  bool side;
593  if( mesh.cells.is_cell_facet_on_surface( c, f, polygon, side ) ) {
594  index_t surface_id { mesh.polygons.surface( polygon ) };
595  side ? out_ << "+" : out_ << "-";
596  out_
597  << region.geomodel().surface( surface_id ).parent( 0 ).name();
598  } else {
599  out_ << "none";
600  }
601  }
602  }
603 
604  void export_model( const GeoModel3D& geomodel )
605  {
606  ringmesh_assert( out_.is_open() );
607  out_ << "MODEL" << EOL;
608  int tface_count = 1;
609 
610  const auto& polygons = geomodel.mesh.polygons;
611  for( auto& cur_interface : geomodel.geol_entities(
612  Interface3D::type_name_static() ) ) {
613  out_ << "SURFACE " << cur_interface.name() << EOL;
614  for( auto s : range( cur_interface.nb_children() ) ) {
615  out_ << "TFACE " << tface_count++ << EOL;
616  auto surface_id = cur_interface.child_gmme( s ).index();
617  out_ << "KEYVERTICES";
618  auto key_polygon_id = polygons.polygon( surface_id, 0 );
619  for( auto v : range( polygons.nb_vertices( key_polygon_id ) ) ) {
620  out_ << " "
621  << vertex_exported_id_[polygons.vertex(
622  { key_polygon_id, v } )];
623  }
624  out_ << EOL;
625  for( auto p : range( polygons.nb_polygons( surface_id ) ) ) {
626  auto polygon_id = polygons.polygon( surface_id, p );
627  out_ << "TRGL";
628  for( auto v : range( polygons.nb_vertices( polygon_id ) ) ) {
629  out_ << " "
630  << vertex_exported_id_[polygons.vertex(
631  { polygon_id, v } )];
632  }
633  out_ << EOL;
634  }
635  }
636  }
637  }
638  void export_model_region( const GeoModel3D& geomodel )
639  {
640  ringmesh_assert( out_.is_open() );
641  for( auto& region : geomodel.regions() ) {
642  out_ << "MODEL_REGION " << region.name() << " ";
643  region.side( 0 ) ? out_ << "+" : out_ << "-";
644  out_ << region.boundary_gmme( 0 ).index() + 1 << EOL;
645  }
646  }
647 
648  void write_no_data_value( index_t nb )
649  {
650  for( auto i : range( nb ) ) {
651  ringmesh_unused( i );
652  out_ << " " << GEO::String::to_string( gocad_no_data_value_ );
653  }
654  }
655 
656  private:
658  std::ofstream out_;
667  std::vector< std::string > numeric_like_vertex_attribute_names_;
668  std::vector< index_t > vertex_attribute_dimensions_;
669  std::vector< std::string > numeric_like_cell_attribute_names_;
670  std::vector< index_t > cell_attribute_dimensions_;
671  std::vector< bool > vertex_exported_;
672  std::vector< bool > atom_exported_;
673  std::vector< index_t > vertex_exported_id_;
674  std::vector< index_t > atom_exported_id_;
675  index_t nb_vertices_exported_ { 1 };
676 
678  const double gocad_no_data_value_ { -99999 };
679  };
680 
681 }
RINGMesh::vec3
vecn< 3 > vec3
Definition: types.h:76
ringmesh_unused
void ringmesh_unused(const T &)
Definition: common.h:105
RINGMesh::region_range
Definition: geomodel_indexing_types.h:157
RINGMesh::GeoModelBuilderTSolid
Builds a meshed GeoModel from a Gocad TSolid (file.so)
Definition: geomodel_builder_gocad.h:545
ringmesh_assert
#define ringmesh_assert(x)
Definition: ringmesh_assert.h:63
RINGMesh::contains
bool contains(const container &in, const T &value)
Definition: algorithm.h:87
RINGMesh
Classes to build GeoModel from various inputs.
Definition: algorithm.h:48
EOL
const char EOL
Definition: io.h:45
ringmesh_assert_not_reached
#define ringmesh_assert_not_reached
Definition: ringmesh_assert.h:64