MOAB: Mesh Oriented datABase  (version 5.2.1)
crop_vol_test.cpp
Go to the documentation of this file.
00001 /**
00002  * This library is free software; you can redistribute it and/or
00003  * modify it under the terms of the GNU Lesser General Public
00004  * License as published by the Free Software Foundation; either
00005  * version 2.1 of the License, or (at your option) any later version.
00006  *
00007  */
00008 /*
00009  * crop_vol_test.cpp
00010  * \brief Driver to create a volume from 2 surfaces (top and bottom), a user defined polyline
00011  * (contour) and a direction
00012  *
00013  */
00014 
00015 #include "moab/Core.hpp"
00016 #include <iostream>
00017 #include <fstream>
00018 #include <string.h>
00019 #include "moab/FBEngine.hpp"
00020 #include "moab/GeomTopoTool.hpp"
00021 #include "TestUtil.hpp"
00022 
00023 std::string filename_top;
00024 std::string filename_bot;
00025 std::string polygon_file_name;
00026 double min_dot = 0.8;
00027 std::string vol_file;
00028 
00029 int number_tests_successful = 0;
00030 int number_tests_failed     = 0;
00031 
00032 using namespace moab;
00033 
00034 ErrorCode volume_test( FBEngine* pFacet );
00035 
00036 void handle_error_code( ErrorCode rv, int& number_failed, int& number_successful )
00037 {
00038     if( rv == MB_SUCCESS )
00039     {
00040         std::cout << "Success";
00041         number_successful++;
00042     }
00043     else
00044     {
00045         std::cout << "Failure";
00046         number_failed++;
00047     }
00048 }
00049 
00050 int main( int argc, char* argv[] )
00051 {
00052     filename_bot      = TestDir + "/BedCrop2.h5m";
00053     filename_top      = TestDir + "/SECrop2.h5m";  //"/polyPB.txt";
00054     polygon_file_name = TestDir + "/poly14.txt";
00055     vol_file          = "volIce.h5m";  // output
00056 
00057     number_tests_successful = 0;
00058     bool remove_output      = true;
00059 
00060     if( argc == 1 )
00061     {
00062         std::cout << "Using default input files: " << filename_bot << " " << filename_top << " " << polygon_file_name
00063                   << " " << min_dot << " " << vol_file << std::endl;
00064         std::cout << "    default output file: " << vol_file << " will be deleted \n";
00065     }
00066     else if( argc == 6 )
00067     {
00068         remove_output     = false;
00069         filename_bot      = argv[1];
00070         filename_top      = argv[2];
00071         polygon_file_name = argv[3];
00072         min_dot           = atof( argv[4] );
00073         vol_file          = argv[5];
00074     }
00075     else
00076     {
00077         std::cout << " wrong number of arguments\n";
00078         return 1;  // fail
00079     }
00080     Core mbcore;
00081     Interface* mb = &mbcore;
00082 
00083     ErrorCode rval = mb->load_file( filename_bot.c_str() );MB_CHK_SET_ERR( rval, "failed to load bed file" );
00084 
00085     rval = mb->load_file( filename_top.c_str() );  // load the second face (we know we have one face in each file)
00086     MB_CHK_SET_ERR( rval, "failed to load top file" );
00087 
00088     FBEngine* pFacet = new FBEngine( mb, NULL, true );  // smooth facetting, no OBB tree passed
00089 
00090     if( !pFacet ) return 1;  // error
00091 
00092     // should the init be part of constructor or not?
00093     // this is where the obb tree is constructed, and smooth faceting initialized, too.
00094     rval = pFacet->Init();MB_CHK_SET_ERR( rval, "failed to initialize smoothing" );
00095 
00096     std::cout << "volume creation test: ";
00097     rval = volume_test( pFacet );
00098     handle_error_code( rval, number_tests_failed, number_tests_successful );
00099     std::cout << "\n";
00100 
00101     if( remove_output ) { remove( vol_file.c_str() ); }
00102     return number_tests_failed;
00103 }
00104 
00105 ErrorCode volume_test( FBEngine* pFacet )
00106 {
00107     // there should be 2 faces in the model
00108     std::ifstream datafile( polygon_file_name.c_str(), std::ifstream::in );
00109     if( !datafile )
00110     {
00111         std::cout << "can't read file\n";
00112         return MB_FAILURE;
00113     }
00114     //
00115     char temp[100];
00116     double direction[3];  // normalized
00117     double gridSize;
00118     datafile.getline( temp, 100 );  // first line
00119 
00120     // get direction and mesh size along polygon segments, from file
00121     sscanf( temp, " %lf %lf %lf %lf ", direction, direction + 1, direction + 2, &gridSize );
00122     // NORMALIZE(direction);// just to be sure
00123 
00124     std::vector< double > xyz;
00125     while( !datafile.eof() )
00126     {
00127         datafile.getline( temp, 100 );
00128         // int id = 0;
00129         double x, y, z;
00130         int nr = sscanf( temp, "%lf %lf %lf", &x, &y, &z );
00131         if( nr == 3 )
00132         {
00133             xyz.push_back( x );
00134             xyz.push_back( y );
00135             xyz.push_back( z );
00136         }
00137     }
00138     int sizePolygon = (int)xyz.size() / 3;
00139     if( sizePolygon < 3 )
00140     {
00141         std::cerr << " Not enough points in the polygon" << std::endl;
00142         return MB_FAILURE;
00143     }
00144     EntityHandle root_set;
00145     ErrorCode rval = pFacet->getRootSet( &root_set );MB_CHK_SET_ERR( rval, "ERROR : getRootSet failed!" );
00146 
00147     int top = 2;  //  iBase_FACE;
00148 
00149     Range faces;
00150     rval = pFacet->getEntities( root_set, top, faces );MB_CHK_SET_ERR( rval, "Failed to get faces in volume_test." );
00151 
00152     if( 2 != faces.size() )
00153     {
00154         std::cerr << "expected 2 faces, top and bottom\n";
00155         return MB_FAILURE;
00156     }
00157     EntityHandle newFace1;  // first test is with closed surface
00158     rval = pFacet->split_surface_with_direction( faces[0], xyz, direction, /*closed*/ 1, min_dot, newFace1 );MB_CHK_SET_ERR( rval, "Failed to crop first face." );
00159 
00160     EntityHandle newFace2;  // first test is with closed surface
00161     rval = pFacet->split_surface_with_direction( faces[1], xyz, direction, /*closed*/ 1, min_dot, newFace2 );MB_CHK_SET_ERR( rval, "Failed to crop second face." );
00162 
00163     // here we make the assumption that the edges, vertices are matching fine...
00164     EntityHandle volume;
00165     rval = pFacet->create_volume_with_direction( newFace1, newFace2, direction, volume );MB_CHK_SET_ERR( rval, "Failed to create volume." );
00166 
00167     Interface* mb = pFacet->moab_instance();
00168     pFacet->delete_smooth_tags();
00169     // duplicate -- extract a new geom topo tool
00170     GeomTopoTool* gtt       = pFacet->get_gtt();
00171     GeomTopoTool* duplicate = NULL;
00172     std::vector< EntityHandle > gents;
00173     gents.push_back( volume );
00174     rval = gtt->duplicate_model( duplicate, &gents );MB_CHK_SET_ERR( rval, "Failed to extract volume." );
00175     EntityHandle newRootSet = duplicate->get_root_model_set();
00176     delete pFacet;
00177     pFacet = NULL;  // try not to write the obb tree
00178     delete duplicate;
00179     rval = mb->write_file( vol_file.c_str(), NULL, NULL, &newRootSet, 1 );
00180 
00181     return rval;
00182 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines