Mesh Oriented datABase  (version 5.4.1)
Array-based unstructured mesh datastructure
ReadCGM.cpp
Go to the documentation of this file.
00001 /**
00002  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
00003  * storing and accessing finite element mesh data.
00004  *
00005  * Copyright 2004 Sandia Corporation.  Under the terms of Contract
00006  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
00007  * retains certain rights in this software.
00008  *
00009  * This library is free software; you can redistribute it and/or
00010  * modify it under the terms of the GNU Lesser General Public
00011  * License as published by the Free Software Foundation; either
00012  * version 2.1 of the License, or (at your option) any later version.
00013  *
00014  */
00015 #ifdef WIN32
00016 #pragma warning( disable : 4786 )
00017 #endif
00018 
00019 #include "CGMConfig.h"
00020 #include "GeometryQueryTool.hpp"
00021 #include "ModelQueryEngine.hpp"
00022 #include "RefEntityName.hpp"
00023 #include "GMem.hpp"
00024 
00025 #include "RefGroup.hpp"
00026 #include "RefVolume.hpp"
00027 #include "RefFace.hpp"
00028 #include "RefEdge.hpp"
00029 #include "RefVertex.hpp"
00030 
00031 #include "SenseEntity.hpp"
00032 #include "Surface.hpp"
00033 #include "Curve.hpp"
00034 #include "Body.hpp"
00035 #include "InitCGMA.hpp"
00036 
00037 #include "moab/Core.hpp"
00038 #include "moab/Interface.hpp"
00039 #include "moab/Range.hpp"
00040 #include "MBTagConventions.hpp"
00041 #include "moab/FileOptions.hpp"
00042 
00043 #include "moab/GeomTopoTool.hpp"
00044 
00045 #include "CubitCompat.hpp"
00046 
00047 #include <cstdio>
00048 #include <algorithm>
00049 #include <cassert>
00050 #include <iostream>
00051 
00052 #include "ReadCGM.hpp"
00053 #include "moab/ReadUtilIface.hpp"
00054 
00055 namespace moab
00056 {
00057 
00058 #define GF_CUBIT_FILE_TYPE    "CUBIT"
00059 #define GF_STEP_FILE_TYPE     "STEP"
00060 #define GF_IGES_FILE_TYPE     "IGES"
00061 #define GF_OCC_BREP_FILE_TYPE "OCC"
00062 #define GF_FACET_FILE_TYPE    "FACET"
00063 
00064 ReaderIface* ReadCGM::factory( Interface* iface )
00065 {
00066     return new ReadCGM( iface );
00067 }
00068 
00069 ReadCGM::ReadCGM( Interface* impl )
00070     : geom_tag( 0 ), id_tag( 0 ), name_tag( 0 ), category_tag( 0 ), faceting_tol_tag( 0 ), geometry_resabs_tag( 0 )
00071 {
00072     assert( NULL != impl );
00073     mdbImpl    = impl;
00074     myGeomTool = new GeomTopoTool( impl );
00075     impl->query_interface( readUtilIface );
00076     assert( NULL != readUtilIface );
00077 
00078     // initialise counters
00079     failed_curve_count   = 0;
00080     failed_surface_count = 0;
00081 
00082     ErrorCode rval;
00083 
00084     // get some tag handles
00085     int negone = -1, zero = 0 /*, negonearr[] = {-1, -1, -1, -1}*/;
00086     rval = mdbImpl->tag_get_handle( GEOM_DIMENSION_TAG_NAME, 1, MB_TYPE_INTEGER, geom_tag, MB_TAG_SPARSE | MB_TAG_CREAT,
00087                                     &negone );
00088     assert( !rval );
00089     id_tag = mdbImpl->globalId_tag();
00090 
00091     rval =
00092         mdbImpl->tag_get_handle( NAME_TAG_NAME, NAME_TAG_SIZE, MB_TYPE_OPAQUE, name_tag, MB_TAG_SPARSE | MB_TAG_CREAT );
00093     assert( !rval );
00094 
00095     rval = mdbImpl->tag_get_handle( CATEGORY_TAG_NAME, CATEGORY_TAG_SIZE, MB_TYPE_OPAQUE, category_tag,
00096                                     MB_TAG_SPARSE | MB_TAG_CREAT );
00097     assert( !rval );
00098     rval = mdbImpl->tag_get_handle( "FACETING_TOL", 1, MB_TYPE_DOUBLE, faceting_tol_tag, MB_TAG_SPARSE | MB_TAG_CREAT );
00099     assert( !rval );
00100     rval = mdbImpl->tag_get_handle( "GEOMETRY_RESABS", 1, MB_TYPE_DOUBLE, geometry_resabs_tag,
00101                                     MB_TAG_SPARSE | MB_TAG_CREAT );
00102     assert( !rval );
00103 #ifdef NDEBUG
00104     if( !rval )
00105     {
00106     };  // Line to avoid compiler warning about variable set but not used
00107 #endif
00108 }
00109 
00110 ReadCGM::~ReadCGM()
00111 {
00112     mdbImpl->release_interface( readUtilIface );
00113     delete myGeomTool;
00114 }
00115 
00116 ErrorCode ReadCGM::read_tag_values( const char* /* file_name */,
00117                                     const char* /* tag_name */,
00118                                     const FileOptions& /* opts */,
00119                                     std::vector< int >& /* tag_values_out */,
00120                                     const SubsetList* /* subset_list */ )
00121 {
00122     return MB_NOT_IMPLEMENTED;
00123 }
00124 
00125 // Sets options passed into ReadCGM::load_file
00126 ErrorCode ReadCGM::set_options( const FileOptions& opts,
00127                                 int& norm_tol,
00128                                 double& faceting_tol,
00129                                 double& len_tol,
00130                                 bool& act_att,
00131                                 bool& verbose_warnings,
00132                                 bool& fatal_on_curves )
00133 {
00134     ErrorCode rval;
00135 
00136     // Default Values
00137     int DEFAULT_NORM         = 5;
00138     double DEFAULT_FACET_TOL = 0.001;
00139     double DEFAULT_LEN_TOL   = 0.0;
00140     act_att                  = true;
00141 
00142     // Check for the options.
00143     if( MB_SUCCESS != opts.get_int_option( "FACET_NORMAL_TOLERANCE", norm_tol ) ) norm_tol = DEFAULT_NORM;
00144 
00145     if( MB_SUCCESS != opts.get_real_option( "FACET_DISTANCE_TOLERANCE", faceting_tol ) )
00146         faceting_tol = DEFAULT_FACET_TOL;
00147 
00148     if( MB_SUCCESS != opts.get_real_option( "MAX_FACET_EDGE_LENGTH", len_tol ) ) len_tol = DEFAULT_LEN_TOL;
00149 
00150     if( MB_SUCCESS == opts.get_null_option( "VERBOSE_CGM_WARNINGS" ) ) verbose_warnings = true;
00151 
00152     if( MB_SUCCESS == opts.get_null_option( "FATAL_ON_CURVES" ) ) fatal_on_curves = true;
00153 
00154     const char* name  = "CGM_ATTRIBS";
00155     const char* value = "no";
00156     rval              = opts.match_option( name, value );
00157     if( MB_SUCCESS == rval ) act_att = false;
00158 
00159     return MB_SUCCESS;
00160 }
00161 
00162 ErrorCode ReadCGM::create_entity_sets( std::map< RefEntity*, EntityHandle > ( &entmap )[5] )
00163 {
00164     ErrorCode rval;
00165     const char geom_categories[][CATEGORY_TAG_SIZE] = { "Vertex\0", "Curve\0", "Surface\0", "Volume\0", "Group\0" };
00166     const char* const names[]                       = { "Vertex", "Curve", "Surface", "Volume" };
00167     DLIList< RefEntity* > entlist;
00168 
00169     for( int dim = 0; dim < 4; dim++ )
00170     {
00171         entlist.clean_out();
00172         GeometryQueryTool::instance()->ref_entity_list( names[dim], entlist, true );
00173         entlist.reset();
00174 
00175         for( int i = entlist.size(); i--; )
00176         {
00177             RefEntity* ent = entlist.get_and_step();
00178             EntityHandle handle;
00179             // Create the new meshset
00180             rval = mdbImpl->create_meshset( dim == 1 ? MESHSET_ORDERED : MESHSET_SET, handle );
00181             if( MB_SUCCESS != rval ) return rval;
00182 
00183             // Map the geom reference entity to the corresponding moab meshset
00184             entmap[dim][ent] = handle;
00185 
00186             // Create tags for the new meshset
00187             rval = mdbImpl->tag_set_data( geom_tag, &handle, 1, &dim );
00188             if( MB_SUCCESS != rval ) return rval;
00189 
00190             int id = ent->id();
00191             rval   = mdbImpl->tag_set_data( id_tag, &handle, 1, &id );
00192             if( MB_SUCCESS != rval ) return rval;
00193 
00194             rval = mdbImpl->tag_set_data( category_tag, &handle, 1, &geom_categories[dim] );
00195             if( MB_SUCCESS != rval ) return rval;
00196         }
00197     }
00198 
00199     return MB_SUCCESS;
00200 }
00201 
00202 ErrorCode ReadCGM::create_topology( std::map< RefEntity*, EntityHandle > ( &entitymap )[5] )
00203 {
00204     ErrorCode rval;
00205     DLIList< RefEntity* > entitylist;
00206     std::map< RefEntity*, EntityHandle >::iterator ci;
00207 
00208     for( int dim = 1; dim < 4; ++dim )
00209     {
00210         for( ci = entitymap[dim].begin(); ci != entitymap[dim].end(); ++ci )
00211         {
00212             entitylist.clean_out();
00213             ci->first->get_child_ref_entities( entitylist );
00214 
00215             entitylist.reset();
00216             for( int i = entitylist.size(); i--; )
00217             {
00218                 RefEntity* ent = entitylist.get_and_step();
00219                 EntityHandle h = entitymap[dim - 1][ent];
00220                 rval           = mdbImpl->add_parent_child( ci->second, h );
00221                 if( MB_SUCCESS != rval ) return rval;
00222             }
00223         }
00224     }
00225 
00226     return MB_SUCCESS;
00227 }
00228 
00229 ErrorCode ReadCGM::store_surface_senses( std::map< RefEntity*, EntityHandle >& surface_map,
00230                                          std::map< RefEntity*, EntityHandle >& volume_map )
00231 {
00232     ErrorCode rval;
00233     std::map< RefEntity*, EntityHandle >::iterator ci;
00234 
00235     for( ci = surface_map.begin(); ci != surface_map.end(); ++ci )
00236     {
00237         RefFace* face                = (RefFace*)( ci->first );
00238         BasicTopologyEntity *forward = 0, *reverse = 0;
00239         for( SenseEntity* cf = face->get_first_sense_entity_ptr(); cf; cf = cf->next_on_bte() )
00240         {
00241             BasicTopologyEntity* vol = cf->get_parent_basic_topology_entity_ptr();
00242             // Allocate vol to the proper topology entity (forward or reverse)
00243             if( cf->get_sense() == CUBIT_UNKNOWN || cf->get_sense() != face->get_surface_ptr()->bridge_sense() )
00244             {
00245                 // Check that each surface has a sense for only one volume
00246                 if( reverse )
00247                 {
00248                     std::cout << "Surface " << face->id() << " has reverse sense "
00249                               << "with multiple volume " << reverse->id() << " and "
00250                               << "volume " << vol->id() << std::endl;
00251                     return MB_FAILURE;
00252                 }
00253                 reverse = vol;
00254             }
00255             if( cf->get_sense() == CUBIT_UNKNOWN || cf->get_sense() == face->get_surface_ptr()->bridge_sense() )
00256             {
00257                 // Check that each surface has a sense for only one volume
00258                 if( forward )
00259                 {
00260                     std::cout << "Surface " << face->id() << " has forward sense "
00261                               << "with multiple volume " << forward->id() << " and "
00262                               << "volume " << vol->id() << std::endl;
00263                     return MB_FAILURE;
00264                 }
00265                 forward = vol;
00266             }
00267         }
00268 
00269         if( forward )
00270         {
00271             rval = myGeomTool->set_sense( ci->second, volume_map[forward], SENSE_FORWARD );
00272             if( MB_SUCCESS != rval ) return rval;
00273         }
00274         if( reverse )
00275         {
00276             rval = myGeomTool->set_sense( ci->second, volume_map[reverse], SENSE_REVERSE );
00277             if( MB_SUCCESS != rval ) return rval;
00278         }
00279     }
00280 
00281     return MB_SUCCESS;
00282 }
00283 
00284 ErrorCode ReadCGM::store_curve_senses( std::map< RefEntity*, EntityHandle >& curve_map,
00285                                        std::map< RefEntity*, EntityHandle >& surface_map )
00286 {
00287     ErrorCode rval;
00288     std::vector< EntityHandle > ents;
00289     std::vector< int > senses;
00290     std::map< RefEntity*, EntityHandle >::iterator ci;
00291     for( ci = curve_map.begin(); ci != curve_map.end(); ++ci )
00292     {
00293         RefEdge* edge = (RefEdge*)( ci->first );
00294         ents.clear();
00295         senses.clear();
00296         for( SenseEntity* ce = edge->get_first_sense_entity_ptr(); ce; ce = ce->next_on_bte() )
00297         {
00298             BasicTopologyEntity* fac = ce->get_parent_basic_topology_entity_ptr();
00299             EntityHandle face        = surface_map[fac];
00300             if( ce->get_sense() == CUBIT_UNKNOWN || ce->get_sense() != edge->get_curve_ptr()->bridge_sense() )
00301             {
00302                 ents.push_back( face );
00303                 senses.push_back( SENSE_REVERSE );
00304             }
00305             if( ce->get_sense() == CUBIT_UNKNOWN || ce->get_sense() == edge->get_curve_ptr()->bridge_sense() )
00306             {
00307                 ents.push_back( face );
00308                 senses.push_back( SENSE_FORWARD );
00309             }
00310         }
00311 
00312         rval = myGeomTool->set_senses( ci->second, ents, senses );
00313         if( MB_SUCCESS != rval ) return rval;
00314     }
00315     return MB_SUCCESS;
00316 }
00317 
00318 ErrorCode ReadCGM::store_groups( std::map< RefEntity*, EntityHandle > ( &entitymap )[5] )
00319 {
00320     ErrorCode rval;
00321 
00322     // Create entity sets for all ref groups
00323     rval = create_group_entsets( entitymap[4] );
00324     if( rval != MB_SUCCESS ) return rval;
00325 
00326     // Store group names and entities in the mesh
00327     rval = store_group_content( entitymap );
00328     if( rval != MB_SUCCESS ) return rval;
00329 
00330     return MB_SUCCESS;
00331 }
00332 
00333 ErrorCode ReadCGM::create_group_entsets( std::map< RefEntity*, EntityHandle >& group_map )
00334 {
00335     ErrorCode rval;
00336     const char geom_categories[][CATEGORY_TAG_SIZE] = { "Vertex\0", "Curve\0", "Surface\0", "Volume\0", "Group\0" };
00337     DLIList< RefEntity* > entitylist;
00338     // Create entity sets for all ref groups
00339     std::vector< Tag > extra_name_tags;
00340 #if CGM_MAJOR_VERSION > 13
00341     DLIList< CubitString > name_list;
00342 #else
00343     DLIList< CubitString* > name_list;
00344 #endif
00345     entitylist.clean_out();
00346     // Get all entity groups from the CGM model
00347     GeometryQueryTool::instance()->ref_entity_list( "group", entitylist );
00348     entitylist.reset();
00349     // Loop over all groups
00350     for( int i = entitylist.size(); i--; )
00351     {
00352         // Take the next group
00353         RefEntity* grp = entitylist.get_and_step();
00354         name_list.clean_out();
00355 // Get the names of all entities in this group from the solid model
00356 #if CGM_MAJOR_VERSION > 13
00357         RefEntityName::instance()->get_refentity_name( grp, name_list );
00358 #else
00359         // True argument is optional, but for large multi-names situation, it should save
00360         // some cpu time
00361         RefEntityName::instance()->get_refentity_name( grp, name_list, true );
00362 #endif
00363         if( name_list.size() == 0 ) continue;
00364         // Set pointer to first name of the group and set the first name to name1
00365         name_list.reset();
00366 #if CGM_MAJOR_VERSION > 13
00367         CubitString name1 = name_list.get();
00368 #else
00369         CubitString name1 = *name_list.get();
00370 #endif
00371         // Create entity handle for the group
00372         EntityHandle h;
00373         rval = mdbImpl->create_meshset( MESHSET_SET, h );
00374         if( MB_SUCCESS != rval ) return rval;
00375         // Set tag data for the group
00376         char namebuf[NAME_TAG_SIZE];
00377         memset( namebuf, '\0', NAME_TAG_SIZE );
00378         strncpy( namebuf, name1.c_str(), NAME_TAG_SIZE - 1 );
00379         if( name1.length() >= (unsigned)NAME_TAG_SIZE )
00380             std::cout << "WARNING: group name '" << name1.c_str() << "' truncated to '" << namebuf << "'" << std::endl;
00381         rval = mdbImpl->tag_set_data( name_tag, &h, 1, namebuf );
00382         if( MB_SUCCESS != rval ) return MB_FAILURE;
00383 
00384         int id = grp->id();
00385         rval   = mdbImpl->tag_set_data( id_tag, &h, 1, &id );
00386         if( MB_SUCCESS != rval ) return MB_FAILURE;
00387 
00388         rval = mdbImpl->tag_set_data( category_tag, &h, 1, &geom_categories[4] );
00389         if( MB_SUCCESS != rval ) return MB_FAILURE;
00390         // Check for extra group names
00391         if( name_list.size() > 1 )
00392         {
00393             for( int j = extra_name_tags.size(); j < name_list.size(); ++j )
00394             {
00395                 sprintf( namebuf, "EXTRA_%s%d", NAME_TAG_NAME, j );
00396                 Tag t;
00397                 rval =
00398                     mdbImpl->tag_get_handle( namebuf, NAME_TAG_SIZE, MB_TYPE_OPAQUE, t, MB_TAG_SPARSE | MB_TAG_CREAT );
00399                 assert( !rval );
00400                 extra_name_tags.push_back( t );
00401             }
00402             // Add extra group names to the group handle
00403             for( int j = 0; j < name_list.size(); ++j )
00404             {
00405 #if CGM_MAJOR_VERSION > 13
00406                 name1 = name_list.get_and_step();
00407 #else
00408                 name1 = *name_list.get_and_step();
00409 #endif
00410                 memset( namebuf, '\0', NAME_TAG_SIZE );
00411                 strncpy( namebuf, name1.c_str(), NAME_TAG_SIZE - 1 );
00412                 if( name1.length() >= (unsigned)NAME_TAG_SIZE )
00413                     std::cout << "WARNING: group name '" << name1.c_str() << "' truncated to '" << namebuf << "'"
00414                               << std::endl;
00415                 rval = mdbImpl->tag_set_data( extra_name_tags[j], &h, 1, namebuf );
00416                 if( MB_SUCCESS != rval ) return MB_FAILURE;
00417             }
00418         }
00419         // Add the group handle
00420         group_map[grp] = h;
00421     }
00422 
00423     return MB_SUCCESS;
00424 }
00425 
00426 ErrorCode ReadCGM::store_group_content( std::map< RefEntity*, EntityHandle > ( &entitymap )[5] )
00427 {
00428     ErrorCode rval;
00429     DLIList< RefEntity* > entlist;
00430     std::map< RefEntity*, EntityHandle >::iterator ci;
00431     // Store contents for each group
00432     entlist.reset();
00433     for( ci = entitymap[4].begin(); ci != entitymap[4].end(); ++ci )
00434     {
00435         RefGroup* grp = (RefGroup*)( ci->first );
00436         entlist.clean_out();
00437         grp->get_child_ref_entities( entlist );
00438 
00439         Range entities;
00440         while( entlist.size() )
00441         {
00442             RefEntity* ent = entlist.pop();
00443             int dim        = ent->dimension();
00444 
00445             if( dim < 0 )
00446             {
00447                 Body* body;
00448                 if( entitymap[4].find( ent ) != entitymap[4].end() )
00449                 {
00450                     // Child is another group; examine its contents
00451                     entities.insert( entitymap[4][ent] );
00452                 }
00453                 else if( ( body = dynamic_cast< Body* >( ent ) ) != NULL )
00454                 {
00455                     // Child is a CGM Body, which presumably comprises some volumes--
00456                     // extract volumes as if they belonged to group.
00457                     DLIList< RefVolume* > vols;
00458                     body->ref_volumes( vols );
00459                     for( int vi = vols.size(); vi--; )
00460                     {
00461                         RefVolume* vol = vols.get_and_step();
00462                         if( entitymap[3].find( vol ) != entitymap[3].end() )
00463                         {
00464                             entities.insert( entitymap[3][vol] );
00465                         }
00466                         else
00467                         {
00468                             std::cerr << "Warning: CGM Body has orphan RefVolume" << std::endl;
00469                         }
00470                     }
00471                 }
00472                 else
00473                 {
00474                     // Otherwise, warn user.
00475                     std::cerr << "Warning: A dim<0 entity is being ignored by ReadCGM." << std::endl;
00476                 }
00477             }
00478             else if( dim < 4 )
00479             {
00480                 if( entitymap[dim].find( ent ) != entitymap[dim].end() ) entities.insert( entitymap[dim][ent] );
00481             }
00482         }
00483 
00484         if( !entities.empty() )
00485         {
00486             rval = mdbImpl->add_entities( ci->second, entities );
00487             if( MB_SUCCESS != rval ) return MB_FAILURE;
00488         }
00489     }
00490 
00491     return MB_SUCCESS;
00492 }
00493 
00494 void ReadCGM::set_cgm_attributes( bool const act_attributes, bool const verbose )
00495 {
00496     if( act_attributes )
00497     {
00498         CGMApp::instance()->attrib_manager()->set_all_auto_read_flags( act_attributes );
00499         CGMApp::instance()->attrib_manager()->set_all_auto_actuate_flags( act_attributes );
00500     }
00501 
00502     if( !verbose )
00503     {
00504         CGMApp::instance()->attrib_manager()->silent_flag( true );
00505     }
00506 }
00507 
00508 ErrorCode ReadCGM::create_vertices( std::map< RefEntity*, EntityHandle >& vertex_map )
00509 {
00510     ErrorCode rval;
00511     std::map< RefEntity*, EntityHandle >::iterator ci;
00512     for( ci = vertex_map.begin(); ci != vertex_map.end(); ++ci )
00513     {
00514         CubitVector pos  = dynamic_cast< RefVertex* >( ci->first )->coordinates();
00515         double coords[3] = { pos.x(), pos.y(), pos.z() };
00516         EntityHandle vh;
00517         rval = mdbImpl->create_vertex( coords, vh );
00518         if( MB_SUCCESS != rval ) return MB_FAILURE;
00519 
00520         // Add the vertex to its tagged meshset
00521         rval = mdbImpl->add_entities( ci->second, &vh, 1 );
00522         if( MB_SUCCESS != rval ) return MB_FAILURE;
00523 
00524         // Replace the meshset handle with the vertex handle
00525         // This makes adding the vertex to higher dim sets easier
00526         ci->second = vh;
00527     }
00528 
00529     return MB_SUCCESS;
00530 }
00531 
00532 ErrorCode ReadCGM::create_curve_facets( std::map< RefEntity*, EntityHandle >& curve_map,
00533                                         std::map< RefEntity*, EntityHandle >& vertex_map,
00534 #if CGM_MAJOR_VERSION > 12
00535                                         int norm_tol,
00536 #else
00537                                         int /* norm_tol */,
00538 #endif
00539                                         double faceting_tol,
00540                                         bool verbose_warn,
00541                                         bool fatal_on_curves )
00542 {
00543     ErrorCode rval;
00544     CubitStatus s;
00545     // Maximum allowable curve-endpoint proximity warnings
00546     // If this integer becomes negative, then abs(curve_warnings) is the
00547     // number of warnings that were suppressed.
00548     int curve_warnings = 0;
00549 
00550     // Map iterator
00551     std::map< RefEntity*, EntityHandle >::iterator ci;
00552 
00553     // Create geometry for all curves
00554     GMem data;
00555     for( ci = curve_map.begin(); ci != curve_map.end(); ++ci )
00556     {
00557         // Get the start and end points of the curve in the form of a reference edge
00558         RefEdge* edge = dynamic_cast< RefEdge* >( ci->first );
00559         // Get the edge's curve information
00560         Curve* curve = edge->get_curve_ptr();
00561         // Clean out previous curve information
00562         data.clean_out();
00563         // Facet curve according to parameters and CGM version
00564 #if CGM_MAJOR_VERSION > 12
00565         s = edge->get_graphics( data, norm_tol, faceting_tol );
00566 #else
00567         s = edge->get_graphics( data, faceting_tol );
00568 #endif
00569 
00570         if( s != CUBIT_SUCCESS )
00571         {
00572             // if we fatal on curves
00573             if( fatal_on_curves )
00574             {
00575                 std::cout << "Failed to facet the curve " << edge->id() << std::endl;
00576                 return MB_FAILURE;
00577             }
00578             // otherwise record them
00579             else
00580             {
00581                 failed_curve_count++;
00582                 failed_curves.push_back( edge->id() );
00583             }
00584             continue;
00585         }
00586 
00587         std::vector< CubitVector > points;
00588         for( int i = 0; i < data.pointListCount; ++i )
00589             // Add Cubit vertext points to a list
00590             points.push_back( CubitVector( data.point_list()[i].x, data.point_list()[i].y, data.point_list()[i].z ) );
00591 
00592         // Need to reverse data?
00593         if( curve->bridge_sense() == CUBIT_REVERSED ) std::reverse( points.begin(), points.end() );
00594 
00595         // Check for closed curve
00596         RefVertex *start_vtx, *end_vtx;
00597         start_vtx = edge->start_vertex();
00598         end_vtx   = edge->end_vertex();
00599 
00600         // Special case for point curve
00601         if( points.size() < 2 )
00602         {
00603             if( start_vtx != end_vtx || curve->measure() > GEOMETRY_RESABS )
00604             {
00605                 std::cerr << "Warning: No facetting for curve " << edge->id() << std::endl;
00606                 continue;
00607             }
00608             EntityHandle h = vertex_map[start_vtx];
00609             rval           = mdbImpl->add_entities( ci->second, &h, 1 );
00610             if( MB_SUCCESS != rval ) return MB_FAILURE;
00611             continue;
00612         }
00613         // Check to see if the first and last interior vertices are considered to be
00614         // coincident by CUBIT
00615         const bool closed = ( points.front() - points.back() ).length() < GEOMETRY_RESABS;
00616         if( closed != ( start_vtx == end_vtx ) )
00617         {
00618             std::cerr << "Warning: topology and geometry inconsistant for possibly closed curve " << edge->id()
00619                       << std::endl;
00620         }
00621 
00622         // Check proximity of vertices to end coordinates
00623         if( ( start_vtx->coordinates() - points.front() ).length() > GEOMETRY_RESABS ||
00624             ( end_vtx->coordinates() - points.back() ).length() > GEOMETRY_RESABS )
00625         {
00626 
00627             curve_warnings--;
00628             if( curve_warnings >= 0 || verbose_warn )
00629             {
00630                 std::cerr << "Warning: vertices not at ends of curve " << edge->id() << std::endl;
00631                 if( curve_warnings == 0 && !verbose_warn )
00632                 {
00633                     std::cerr << "         further instances of this warning will be suppressed..." << std::endl;
00634                 }
00635             }
00636         }
00637         // Create interior points
00638         std::vector< EntityHandle > verts, edges;
00639         verts.push_back( vertex_map[start_vtx] );
00640         for( size_t i = 1; i < points.size() - 1; ++i )
00641         {
00642             double coords[] = { points[i].x(), points[i].y(), points[i].z() };
00643             EntityHandle h;
00644             // Create vertex entity
00645             rval = mdbImpl->create_vertex( coords, h );
00646             if( MB_SUCCESS != rval ) return MB_FAILURE;
00647             verts.push_back( h );
00648         }
00649         verts.push_back( vertex_map[end_vtx] );
00650 
00651         // Create edges
00652         for( size_t i = 0; i < verts.size() - 1; ++i )
00653         {
00654             EntityHandle h;
00655             rval = mdbImpl->create_element( MBEDGE, &verts[i], 2, h );
00656             if( MB_SUCCESS != rval ) return MB_FAILURE;
00657             edges.push_back( h );
00658         }
00659 
00660         // If closed, remove duplicate
00661         if( verts.front() == verts.back() ) verts.pop_back();
00662         // Add entities to the curve meshset from entitymap
00663         rval = mdbImpl->add_entities( ci->second, &verts[0], verts.size() );
00664         if( MB_SUCCESS != rval ) return MB_FAILURE;
00665         rval = mdbImpl->add_entities( ci->second, &edges[0], edges.size() );
00666         if( MB_SUCCESS != rval ) return MB_FAILURE;
00667     }
00668 
00669     if( !verbose_warn && curve_warnings < 0 )
00670     {
00671         std::cerr << "Suppressed " << -curve_warnings << " 'vertices not at ends of curve' warnings." << std::endl;
00672         std::cerr << "To see all warnings, use reader param VERBOSE_CGM_WARNINGS." << std::endl;
00673     }
00674 
00675     return MB_SUCCESS;
00676 }
00677 
00678 ErrorCode ReadCGM::create_surface_facets( std::map< RefEntity*, EntityHandle >& surface_map,
00679                                           std::map< RefEntity*, EntityHandle >& vertex_map,
00680                                           int norm_tol,
00681                                           double facet_tol,
00682                                           double length_tol )
00683 {
00684     ErrorCode rval;
00685     std::map< RefEntity*, EntityHandle >::iterator ci;
00686     CubitStatus s;
00687 #if( ( CGM_MAJOR_VERSION == 14 && CGM_MINOR_VERSION > 2 ) || CGM_MAJOR_VERSION >= 15 )
00688     DLIList< TopologyEntity* > me_list;
00689 #else
00690     DLIList< ModelEntity* > me_list;
00691 #endif
00692 
00693     GMem data;
00694     // Create geometry for all surfaces
00695     for( ci = surface_map.begin(); ci != surface_map.end(); ++ci )
00696     {
00697         RefFace* face = dynamic_cast< RefFace* >( ci->first );
00698 
00699         data.clean_out();
00700         s = face->get_graphics( data, norm_tol, facet_tol, length_tol );
00701 
00702         if( CUBIT_SUCCESS != s ) return MB_FAILURE;
00703 
00704         // Declare array of all vertex handles
00705         std::vector< EntityHandle > verts( data.pointListCount, 0 );
00706 
00707         // Get list of geometric vertices in surface
00708         me_list.clean_out();
00709         ModelQueryEngine::instance()->query_model( *face, DagType::ref_vertex_type(), me_list );
00710 
00711         // For each geometric vertex, find a single coincident point in facets
00712         // Otherwise, print a warning
00713         for( int i = me_list.size(); i--; )
00714         {
00715             // Assign geometric vertex
00716             RefVertex* vtx  = dynamic_cast< RefVertex* >( me_list.get_and_step() );
00717             CubitVector pos = vtx->coordinates();
00718 
00719             for( int j = 0; j < data.pointListCount; ++j )
00720             {
00721                 // Assign facet vertex
00722                 CubitVector vpos( data.point_list()[j].x, data.point_list()[j].y, data.point_list()[j].z );
00723                 // Check to see if they are considered coincident
00724                 if( ( pos - vpos ).length_squared() < GEOMETRY_RESABS * GEOMETRY_RESABS )
00725                 {
00726                     // If this facet vertex has already been found coincident, print warning
00727                     if( verts[j] ) std::cerr << "Warning: Coincident vertices in surface " << face->id() << std::endl;
00728                     // If a coincidence is found, keep track of it in the verts vector
00729                     verts[j] = vertex_map[vtx];
00730                     break;
00731                 }
00732             }
00733         }
00734 
00735         // Now create vertices for the remaining points in the facetting
00736         for( int i = 0; i < data.pointListCount; ++i )
00737         {
00738             if( verts[i] )  // If a geometric vertex
00739                 continue;
00740             double coords[] = { data.point_list()[i].x, data.point_list()[i].y, data.point_list()[i].z };
00741             // Return vertex handle to verts to fill in all remaining facet
00742             // vertices
00743             rval = mdbImpl->create_vertex( coords, verts[i] );
00744             if( MB_SUCCESS != rval ) return rval;
00745         }
00746 
00747         // record the failures for information
00748         if( data.fListCount == 0 )
00749         {
00750             failed_surface_count++;
00751             failed_surfaces.push_back( face->id() );
00752         }
00753 
00754         // Now create facets
00755         Range facets;
00756         std::vector< EntityHandle > corners;
00757         for( int i = 0; i < data.fListCount; i += data.facet_list()[i] + 1 )
00758         {
00759             // Get number of facet verts
00760             int* facet = data.facet_list() + i;
00761             corners.resize( *facet );
00762             for( int j = 1; j <= *facet; ++j )
00763             {
00764                 if( facet[j] >= (int)verts.size() )
00765                 {
00766                     std::cerr << "ERROR: Invalid facet data for surface " << face->id() << std::endl;
00767                     return MB_FAILURE;
00768                 }
00769                 corners[j - 1] = verts[facet[j]];
00770             }
00771             EntityType type;
00772             if( *facet == 3 )
00773                 type = MBTRI;
00774             else
00775             {
00776                 std::cerr << "Warning: non-triangle facet in surface " << face->id() << std::endl;
00777                 std::cerr << "  entity has " << *facet << " edges" << std::endl;
00778                 if( *facet == 4 )
00779                     type = MBQUAD;
00780                 else
00781                     type = MBPOLYGON;
00782             }
00783 
00784             // if (surf->bridge_sense() == CUBIT_REVERSED)
00785             // std::reverse(corners.begin(), corners.end());
00786 
00787             EntityHandle h;
00788             rval = mdbImpl->create_element( type, &corners[0], corners.size(), h );
00789             if( MB_SUCCESS != rval ) return MB_FAILURE;
00790 
00791             facets.insert( h );
00792         }
00793 
00794         // Add vertices and facets to surface set
00795         rval = mdbImpl->add_entities( ci->second, &verts[0], verts.size() );
00796         if( MB_SUCCESS != rval ) return MB_FAILURE;
00797         rval = mdbImpl->add_entities( ci->second, facets );
00798         if( MB_SUCCESS != rval ) return MB_FAILURE;
00799     }
00800 
00801     return MB_SUCCESS;
00802 }
00803 
00804 // Copy geometry into mesh database
00805 ErrorCode ReadCGM::load_file( const char* cgm_file_name,
00806                               const EntityHandle* file_set,
00807                               const FileOptions& opts,
00808                               const ReaderIface::SubsetList* subset_list,
00809                               const Tag* /*file_id_tag*/ )
00810 {
00811     // Blocks_to_load and num_blocks are ignored.
00812     ErrorCode rval;
00813 
00814     if( subset_list )
00815     {
00816         MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "Reading subset of files not supported for CGM data" );
00817     }
00818 
00819     int norm_tol;
00820     double faceting_tol;
00821     double len_tol;
00822     bool act_att          = true;
00823     bool verbose_warnings = false;
00824     bool fatal_on_curves  = false;
00825 
00826     rval = set_options( opts, norm_tol, faceting_tol, len_tol, act_att, verbose_warnings, fatal_on_curves );
00827     if( MB_SUCCESS != rval ) return rval;
00828 
00829     // Always tag with the faceting_tol and geometry absolute resolution
00830     // If file_set is defined, use that, otherwise (file_set == NULL) tag the interface
00831     EntityHandle set = file_set ? *file_set : 0;
00832     rval             = mdbImpl->tag_set_data( faceting_tol_tag, &set, 1, &faceting_tol );
00833     if( MB_SUCCESS != rval ) return rval;
00834 
00835     rval = mdbImpl->tag_set_data( geometry_resabs_tag, &set, 1, &GEOMETRY_RESABS );
00836     if( MB_SUCCESS != rval ) return rval;
00837 
00838     // Initialize CGM
00839     InitCGMA::initialize_cgma();
00840 
00841     // Determine CGM settings and amount of output
00842     set_cgm_attributes( act_att, verbose_warnings );
00843 
00844     CubitStatus s;
00845 
00846     // Get CGM file type
00847     const char* file_type = 0;
00848     file_type             = get_geom_file_type( cgm_file_name );
00849     if( !file_type || !strcmp( file_type, "CUBIT" ) ) return MB_FAILURE;
00850 
00851     s = CubitCompat_import_solid_model( cgm_file_name, file_type );
00852     if( CUBIT_SUCCESS != s )
00853     {
00854         MB_SET_ERR( MB_FAILURE, cgm_file_name << ": Failed to read file of type \"" << file_type << "\"" );
00855     }
00856 
00857     // Create entity sets for all geometric entities
00858     std::map< RefEntity*, EntityHandle > entmap[5];  // One for each dim, and one for groups
00859 
00860     rval = create_entity_sets( entmap );
00861     if( rval != MB_SUCCESS ) return rval;
00862 
00863     // Create topology for all geometric entities
00864     rval = create_topology( entmap );
00865     if( rval != MB_SUCCESS ) return rval;
00866 
00867     // Store CoFace senses
00868     rval = store_surface_senses( entmap[2], entmap[3] );
00869     if( rval != MB_SUCCESS ) return rval;
00870 
00871     // Store CoEdge senses
00872     rval = store_curve_senses( entmap[1], entmap[2] );
00873     if( rval != MB_SUCCESS ) return rval;
00874 
00875     // Get group information and store it in the mesh
00876     rval = store_groups( entmap );
00877     if( rval != MB_SUCCESS ) return rval;
00878 
00879     // Done with volumes and groups
00880     entmap[3].clear();
00881     entmap[4].clear();
00882 
00883     // Create geometry for all vertices and replace
00884     rval = create_vertices( entmap[0] );
00885     if( rval != MB_SUCCESS ) return rval;
00886 
00887     // Create facets for all curves
00888     rval = create_curve_facets( entmap[1], entmap[0], norm_tol, faceting_tol, verbose_warnings, fatal_on_curves );
00889     if( rval != MB_SUCCESS ) return rval;
00890 
00891     // Create facets for surfaces
00892     rval = create_surface_facets( entmap[2], entmap[0], norm_tol, faceting_tol, len_tol );
00893     if( rval != MB_SUCCESS ) return rval;
00894 
00895     // print the fail information
00896     dump_fail_counts();
00897 
00898     return MB_SUCCESS;
00899 }
00900 
00901 // return the number of curves that failed to facet
00902 int ReadCGM::get_failed_curve_count()
00903 {
00904     return failed_curve_count;
00905 }
00906 
00907 // return the number of surfaces that failed to facet
00908 int ReadCGM::get_failed_surface_count()
00909 {
00910     return failed_surface_count;
00911 }
00912 
00913 void ReadCGM::dump_fail_counts()
00914 {
00915     std::cout << "***** Faceting Summary Information *****" << std::endl;
00916     std::cout << "----- Curve Fail Information -----" << std::endl;
00917     std::cout << "There were " << failed_curve_count << " curves that could not be faceted." << std::endl;
00918 
00919     if( failed_curve_count > 0 )
00920     {
00921         std::cout << "The curves were ";
00922         for( int i = 0; i < failed_curve_count; i++ )
00923         {
00924             std::cout << failed_curves[i] << " ";
00925             if( ( i % 10 == 0 ) & ( i > 0 ) ) std::cout << std::endl;
00926         }
00927     }
00928     std::cout << std::endl;
00929     std::cout << "----- Facet Fail Information -----" << std::endl;
00930     std::cout << "There were " << failed_surface_count << " surfaces that could not be faceted." << std::endl;
00931     if( failed_surface_count > 0 )
00932     {
00933         std::cout << "The surfaces were ";
00934         for( int i = 0; i < failed_surface_count; i++ )
00935         {
00936             std::cout << failed_surfaces[i] << " ";
00937             if( ( i % 10 == 0 ) & ( i > 0 ) ) std::cout << std::endl;
00938         }
00939     }
00940     std::cout << std::endl;
00941     std::cout << "***** End of Faceting Summary Information *****" << std::endl;
00942     return;
00943 }
00944 
00945 const char* ReadCGM::get_geom_file_type( const char* name )
00946 {
00947     FILE* file;
00948     const char* result = 0;
00949 
00950     file = fopen( name, "r" );
00951     if( file )
00952     {
00953         result = get_geom_fptr_type( file );
00954         fclose( file );
00955     }
00956 
00957     return result;
00958 }
00959 
00960 const char* ReadCGM::get_geom_fptr_type( FILE* file )
00961 {
00962     static const char* CUBIT_NAME = GF_CUBIT_FILE_TYPE;
00963     static const char* STEP_NAME  = GF_STEP_FILE_TYPE;
00964     static const char* IGES_NAME  = GF_IGES_FILE_TYPE;
00965     static const char* BREP_NAME  = GF_OCC_BREP_FILE_TYPE;
00966     static const char* FACET_NAME = GF_FACET_FILE_TYPE;
00967 
00968     if( is_cubit_file( file ) )
00969         return CUBIT_NAME;
00970     else if( is_step_file( file ) )
00971         return STEP_NAME;
00972     else if( is_iges_file( file ) )
00973         return IGES_NAME;
00974     else if( is_occ_brep_file( file ) )
00975         return BREP_NAME;
00976     else if( is_facet_file( file ) )
00977         return FACET_NAME;
00978     else
00979         return NULL;
00980 }
00981 
00982 int ReadCGM::is_cubit_file( FILE* file )
00983 {
00984     unsigned char buffer[4];
00985     return !fseek( file, 0, SEEK_SET ) && fread( buffer, 4, 1, file ) && !memcmp( buffer, "CUBE", 4 );
00986 }
00987 
00988 int ReadCGM::is_step_file( FILE* file )
00989 {
00990     unsigned char buffer[9];
00991     return !fseek( file, 0, SEEK_SET ) && fread( buffer, 9, 1, file ) && !memcmp( buffer, "ISO-10303", 9 );
00992 }
00993 
00994 int ReadCGM::is_iges_file( FILE* file )
00995 {
00996     unsigned char buffer[10];
00997     return !fseek( file, 72, SEEK_SET ) && fread( buffer, 10, 1, file ) && !memcmp( buffer, "S      1", 8 );
00998 }
00999 
01000 int ReadCGM::is_occ_brep_file( FILE* file )
01001 {
01002     unsigned char buffer[6];
01003     return !fseek( file, 0, SEEK_SET ) && fread( buffer, 6, 1, file ) && !memcmp( buffer, "DBRep_", 6 );
01004 }
01005 int ReadCGM::is_facet_file( FILE* file )
01006 {
01007     unsigned char buffer[10];
01008     return !fseek( file, 0, SEEK_SET ) && fread( buffer, 10, 1, file ) && !memcmp( buffer, "MESH_BASED", 10 );
01009 }
01010 
01011 void ReadCGM::tokenize( const std::string& str, std::vector< std::string >& tokens, const char* delimiters )
01012 {
01013     std::string::size_type last = str.find_first_not_of( delimiters, 0 );
01014     std::string::size_type pos  = str.find_first_of( delimiters, last );
01015     while( std::string::npos != pos && std::string::npos != last )
01016     {
01017         tokens.push_back( str.substr( last, pos - last ) );
01018         last = str.find_first_not_of( delimiters, pos );
01019         pos  = str.find_first_of( delimiters, last );
01020         if( std::string::npos == pos ) pos = str.size();
01021     }
01022 }
01023 
01024 }  // namespace moab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines