moab
|
00001 /* 00002 * intx_in_plane_test.cpp 00003 * 00004 * Created on: Oct 4, 2012 00005 */ 00006 #include <iostream> 00007 #include <sstream> 00008 #include <time.h> 00009 #include <stdlib.h> 00010 #include <stdio.h> 00011 #include <string.h> 00012 #include "moab/Core.hpp" 00013 #include "moab/Interface.hpp" 00014 #include "Intx2MeshInPlane.hpp" 00015 #include <math.h> 00016 00017 #define STRINGIFY_(X) #X 00018 #define STRINGIFY(X) STRINGIFY_(X) 00019 00020 using namespace moab; 00021 00022 int main(int argc, char* argv[]) 00023 { 00024 // check command line arg 00025 const char *filename_mesh1 = STRINGIFY(SRCDIR) "/m1.vtk"; 00026 const char *filename_mesh2 = STRINGIFY(SRCDIR) "/m2.vtk"; 00027 const char *newFile = "intx1.vtk"; 00028 const char *edgesFile = "polyWithEdges.vtk"; 00029 if (argc == 4) 00030 { 00031 filename_mesh1 = argv[1]; 00032 filename_mesh2 = argv[2]; 00033 newFile = argv[3]; 00034 } 00035 else 00036 { 00037 printf("Usage: %s <mesh_filename1> <mesh_filename2> <newFile>\n", argv[0]); 00038 if (argc != 1) 00039 return 1; 00040 printf("No files specified. Defaulting to: %s %s %s\n", 00041 filename_mesh1, filename_mesh2, newFile); 00042 } 00043 00044 // read meshes in 2 file sets 00045 ErrorCode rval = MB_SUCCESS; 00046 Core moab; 00047 Interface* mb = &moab;// global 00048 EntityHandle sf1, sf2; 00049 rval = mb->create_meshset(MESHSET_SET, sf1); 00050 if (MB_SUCCESS != rval) 00051 return 1; 00052 rval = mb->create_meshset(MESHSET_SET, sf2); 00053 if (MB_SUCCESS != rval) 00054 return 1; 00055 rval=mb->load_file(filename_mesh1, &sf1); 00056 if (MB_SUCCESS != rval) 00057 return 1; 00058 rval=mb->load_file(filename_mesh2, &sf2); 00059 if (MB_SUCCESS != rval) 00060 return 1; 00061 00062 00063 EntityHandle outputSet; 00064 rval = mb->create_meshset(MESHSET_SET, outputSet); 00065 if (MB_SUCCESS != rval) 00066 return 1; 00067 Intx2MeshInPlane worker(mb); 00068 00069 worker.SetErrorTolerance( 1.e-5); 00070 //worker.enable_debug(); 00071 rval = worker.intersect_meshes(sf1, sf2, outputSet); 00072 if (MB_SUCCESS != rval) 00073 return 1; 00074 rval = mb->write_mesh(newFile, &outputSet, 1); 00075 if (MB_SUCCESS != rval) 00076 return 1; 00077 00078 // retrieve polygons and get adjacent edges 00079 Range polygons; 00080 rval = mb->get_entities_by_type(outputSet, MBPOLYGON, polygons); 00081 if (MB_SUCCESS != rval) 00082 return 1; 00083 00084 Range edges; 00085 rval = mb->get_adjacencies(polygons, 1, true, edges , Interface::UNION); 00086 if (MB_SUCCESS != rval) 00087 return 1; 00088 00089 std::cout << "number of edges:" << edges.size() << "\n"; 00090 // add edges to the output set 00091 rval = mb->add_entities(outputSet, edges); 00092 if (MB_SUCCESS != rval) 00093 return 1; 00094 rval = mb->write_mesh(edgesFile, &outputSet, 1); 00095 if (MB_SUCCESS != rval) 00096 return 1; 00097 return 0; 00098 00099 00100 }