Mesh Oriented datABase
(version 5.4.1)
Array-based unstructured mesh datastructure
|
00001 // Brandon Smith 00002 // January 20, 2009 00003 00004 /* 00005 loop over all surface meshsets 00006 get all quads of surface meshset 00007 loop over all quads 00008 make tris from quad 00009 add tris to the surface meshset 00010 remove quad from surface meshset 00011 delete quad 00012 end loop 00013 end loop 00014 */ 00015 00016 #include "quads_to_tris.hpp" 00017 #include "moab/CartVect.hpp" 00018 using namespace moab; 00019 00020 // Generic function to create two tris from a quad. This can be improved later. 00021 ErrorCode make_tris_from_quad( Interface* MBI, 00022 EntityHandle quad, /* input */ 00023 EntityHandle& tri0, /* output */ 00024 EntityHandle& tri1 /* output */ ) 00025 { 00026 00027 // get connectivity (ordered counterclockwise for 2D elements in MOAB) 00028 ErrorCode result; 00029 const EntityHandle* quad_conn; 00030 int n_verts = 0; 00031 result = MBI->get_connectivity( quad, quad_conn, n_verts ); 00032 assert( 4 == n_verts ); 00033 assert( MB_SUCCESS == result ); 00034 00035 // find length of diagonals 00036 std::vector< CartVect > coords( n_verts ); 00037 result = MBI->get_coords( quad_conn, n_verts, coords[0].array() ); 00038 if( MB_SUCCESS != result ) return result; 00039 CartVect diagA = coords[0] - coords[2]; 00040 double lenA_sqr = diagA.length_squared(); 00041 CartVect diagB = coords[1] - coords[3]; 00042 double lenB_sqr = diagB.length_squared(); 00043 00044 // choose the shortest diagonal 00045 EntityHandle tri0_conn[3], tri1_conn[3]; 00046 if( lenA_sqr < lenB_sqr ) 00047 { 00048 tri0_conn[0] = quad_conn[0]; 00049 tri0_conn[1] = quad_conn[1]; 00050 tri0_conn[2] = quad_conn[2]; 00051 tri1_conn[0] = quad_conn[0]; 00052 tri1_conn[1] = quad_conn[2]; 00053 tri1_conn[2] = quad_conn[3]; 00054 } 00055 else 00056 { 00057 tri0_conn[0] = quad_conn[0]; 00058 tri0_conn[1] = quad_conn[1]; 00059 tri0_conn[2] = quad_conn[3]; 00060 tri1_conn[0] = quad_conn[1]; 00061 tri1_conn[1] = quad_conn[2]; 00062 tri1_conn[2] = quad_conn[3]; 00063 } 00064 00065 // make tris from quad 00066 result = MBI->create_element( MBTRI, tri0_conn, 3, tri0 ); 00067 assert( MB_SUCCESS == result ); 00068 result = MBI->create_element( MBTRI, tri1_conn, 3, tri1 ); 00069 assert( MB_SUCCESS == result ); 00070 00071 return MB_SUCCESS; 00072 } 00073 00074 ErrorCode make_tris_from_quads( Interface* MBI, const Range& quads, Range& tris ) 00075 { 00076 tris.clear(); 00077 for( Range::const_iterator i = quads.begin(); i != quads.end(); ++i ) 00078 { 00079 EntityHandle tri0, tri1; 00080 ErrorCode result = make_tris_from_quad( MBI, *i, tri0, tri1 ); 00081 assert( MB_SUCCESS == result ); 00082 if( MB_SUCCESS != result ) return result; 00083 tris.insert( tri0 ); 00084 tris.insert( tri1 ); 00085 } 00086 return MB_SUCCESS; 00087 } 00088 00089 ErrorCode quads_to_tris( Interface* MBI, EntityHandle input_meshset ) 00090 { 00091 00092 // create a geometry tag to find the surfaces with 00093 ErrorCode result; 00094 Tag geom_tag, id_tag; 00095 result = MBI->tag_get_handle( GEOM_DIMENSION_TAG_NAME, 1, MB_TYPE_INTEGER, geom_tag, MB_TAG_DENSE | MB_TAG_CREAT ); 00096 if( MB_SUCCESS != result ) return result; 00097 00098 // create an id tag to find the surface id with 00099 id_tag = MBI->globalId_tag(); 00100 00101 // get all surface meshsets 00102 Range surface_meshsets; 00103 int dim = 2; 00104 void* input_dim[] = { &dim }; 00105 result = MBI->get_entities_by_type_and_tag( input_meshset, MBENTITYSET, &geom_tag, input_dim, 1, surface_meshsets ); 00106 assert( MB_SUCCESS == result ); 00107 std::cout << surface_meshsets.size() << " surfaces found." << std::endl; 00108 00109 // ****************************************************************** 00110 // Loop over every surface meshset and grab each surface's quads. 00111 // ****************************************************************** 00112 for( Range::iterator i = surface_meshsets.begin(); i != surface_meshsets.end(); ++i ) 00113 { 00114 00115 // get the surface id of the surface meshset 00116 int surf_id = 0; 00117 result = MBI->tag_get_data( id_tag, &( *i ), 1, &surf_id ); 00118 assert( MB_SUCCESS == result ); 00119 std::cout << " Surface " << surf_id << " has "; 00120 00121 // get all quads of the surface 00122 Range quads; 00123 result = MBI->get_entities_by_type( *i, MBQUAD, quads ); 00124 assert( MB_SUCCESS == result ); 00125 std::cout << quads.size() << " quads." << std::endl; 00126 00127 // ****************************************************************** 00128 // For each quad, make two triangles then delete the quad. 00129 // ****************************************************************** 00130 for( Range::iterator j = quads.begin(); j != quads.end(); ++j ) 00131 { 00132 00133 // make the tris 00134 EntityHandle tri0 = 0, tri1 = 0; 00135 result = make_tris_from_quad( MBI, *j, tri0, tri1 ); 00136 assert( MB_SUCCESS == result ); 00137 00138 // add all the tris to the same surface meshset as the quads were inside. 00139 result = MBI->add_entities( *i, &tri0, 1 ); 00140 if( MB_SUCCESS != result ) std::cout << "result=" << result << std::endl; 00141 assert( MB_SUCCESS == result ); 00142 result = MBI->add_entities( *i, &tri1, 1 ); 00143 assert( MB_SUCCESS == result ); 00144 00145 // remove the quad from the surface meshset 00146 result = MBI->remove_entities( *i, &( *j ), 1 ); 00147 assert( MB_SUCCESS == result ); 00148 00149 // delete the quad 00150 result = MBI->delete_entities( &( *j ), 1 ); 00151 assert( MB_SUCCESS == result ); 00152 00153 } // end quad loop 00154 } // end surface meshset loop 00155 return MB_SUCCESS; 00156 }