moab
|
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, /* intput */ 00023 EntityHandle &tri0, /* output */ 00024 EntityHandle &tri1 /* output */) { 00025 00026 // get connectivity (ordered counterclockwise for 2D elements in MOAB) 00027 ErrorCode result; 00028 const EntityHandle *quad_conn; 00029 int n_verts=0; 00030 result = MBI->get_connectivity( quad, quad_conn, n_verts ); 00031 assert( 4 == n_verts ); 00032 assert( MB_SUCCESS == result); 00033 00034 // find length of diagonals 00035 std::vector<CartVect> coords(n_verts); 00036 result = MBI->get_coords( quad_conn, n_verts, coords[0].array() ); 00037 if(MB_SUCCESS != result) return result; 00038 CartVect diagA = coords[0] - coords[2]; 00039 double lenA_sqr= diagA.length_squared(); 00040 CartVect diagB = coords[1] - coords[3]; 00041 double lenB_sqr= diagB.length_squared(); 00042 00043 // choose the shortest diagonal 00044 EntityHandle tri0_conn[3], tri1_conn[3]; 00045 if(lenA_sqr < lenB_sqr) { 00046 tri0_conn[0] = quad_conn[0]; 00047 tri0_conn[1] = quad_conn[1]; 00048 tri0_conn[2] = quad_conn[2]; 00049 tri1_conn[0] = quad_conn[0]; 00050 tri1_conn[1] = quad_conn[2]; 00051 tri1_conn[2] = quad_conn[3]; 00052 } else { 00053 tri0_conn[0] = quad_conn[0]; 00054 tri0_conn[1] = quad_conn[1]; 00055 tri0_conn[2] = quad_conn[3]; 00056 tri1_conn[0] = quad_conn[1]; 00057 tri1_conn[1] = quad_conn[2]; 00058 tri1_conn[2] = quad_conn[3]; 00059 } 00060 00061 // make tris from quad 00062 result = MBI->create_element(MBTRI, tri0_conn, 3, tri0); 00063 assert( MB_SUCCESS == result); 00064 result = MBI->create_element(MBTRI, tri1_conn, 3, tri1); 00065 assert( MB_SUCCESS == result); 00066 00067 return MB_SUCCESS; 00068 } 00069 00070 ErrorCode make_tris_from_quads( Interface *MBI, 00071 const Range quads, 00072 Range &tris ) { 00073 tris.clear(); 00074 for(Range::const_iterator i=quads.begin(); i!=quads.end(); ++i) { 00075 EntityHandle tri0, tri1; 00076 ErrorCode result = make_tris_from_quad( MBI, *i, tri0, tri1 ); 00077 assert(MB_SUCCESS == result); 00078 if (MB_SUCCESS != result) return result; 00079 tris.insert( tri0 ); 00080 tris.insert( tri1 ); 00081 } 00082 return MB_SUCCESS; 00083 } 00084 00085 ErrorCode quads_to_tris( Interface *MBI, EntityHandle input_meshset ) { 00086 00087 // create a geometry tag to find the surfaces with 00088 ErrorCode result; 00089 Tag geom_tag, id_tag; 00090 result = MBI->tag_get_handle( GEOM_DIMENSION_TAG_NAME, 1, MB_TYPE_INTEGER, 00091 geom_tag, MB_TAG_DENSE|MB_TAG_CREAT ); 00092 if (MB_SUCCESS != result) 00093 return result; 00094 00095 // create an id tag to find the surface id with 00096 result = MBI->tag_get_handle( GLOBAL_ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag, 00097 MB_TAG_DENSE|MB_TAG_CREAT ); 00098 if (MB_SUCCESS != result) 00099 return result; 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, 00106 input_dim, 1, surface_meshsets); 00107 assert( MB_SUCCESS == result ); 00108 std::cout << surface_meshsets.size() << " surfaces found." << std::endl; 00109 00110 // ****************************************************************** 00111 // Loop over every surface meshset and grab each surface's quads. 00112 // ****************************************************************** 00113 for( Range::iterator i=surface_meshsets.begin(); i!=surface_meshsets.end(); i++ ) { 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 // make the tris 00133 EntityHandle tri0 = 0, tri1 = 0; 00134 result = make_tris_from_quad( MBI, *j, tri0, tri1 ); 00135 assert( MB_SUCCESS == result ); 00136 00137 // add all the tris to the same surface meshset as the quads were inside. 00138 result = MBI->add_entities( *i, &tri0, 1 ); 00139 if ( MB_SUCCESS != result ) std::cout << "result=" << result << std::endl; 00140 assert( MB_SUCCESS == result); 00141 result = MBI->add_entities( *i, &tri1, 1 ); 00142 assert( MB_SUCCESS == result); 00143 00144 // remove the quad from the surface meshset 00145 result = MBI->remove_entities( *i, &(*j), 1 ); 00146 assert( MB_SUCCESS == result); 00147 00148 // delete the quad 00149 result = MBI->delete_entities( &(*j), 1 ); 00150 assert( MB_SUCCESS == result); 00151 00152 } // end quad loop 00153 } // end surface meshset loop 00154 return MB_SUCCESS; 00155 }