MOAB: Mesh Oriented datABase  (version 5.2.1)
mbaddchunk.cpp
Go to the documentation of this file.
00001 /*
00002  * addncdata.cpp
00003  * this tool will take an existing h5m file and add data from some chunk description files
00004  *  generated from e3sm runs;
00005  * will support mainly showing the chunks in ViSit
00006  *
00007  * example of usage:
00008  * ./mbaddchunk -i penta3d.h5m -n chunks_on_proc.txt -o penta3d_ch.h5m
00009  *
00010  *
00011  * Basically, will output a new h5m file (penta3d_ch.h5m), which has 2 extra tags, corresponding to
00012  * the chunks number and processors it sits on
00013  *
00014  *
00015  *  file  penta3d.h5m is obtained from the pentagons file, with a python script
00016  *
00017  */
00018 
00019 #include "moab/ProgOptions.hpp"
00020 #include "moab/Core.hpp"
00021 
00022 #include <math.h>
00023 #include <sstream>
00024 #include <iostream>
00025 #include <iomanip>
00026 #include <fstream>
00027 
00028 using namespace moab;
00029 using namespace std;
00030 
00031 int main( int argc, char* argv[] )
00032 {
00033 
00034     ProgOptions opts;
00035 
00036     std::string inputfile( "penta3d.h5m" ), outfile( "penta3d_ch.h5m" ), chunkfile_name, gsmapfile;
00037 
00038     opts.addOpt< std::string >( "input,i", "input mesh filename", &inputfile );
00039     opts.addOpt< std::string >( "chunkFile,n", "chunk file from cam run", &chunkfile_name );
00040     opts.addOpt< std::string >( "gsMAPfile,g", "gsmap file", &gsmapfile );
00041 
00042     opts.addOpt< std::string >( "output,o", "output mesh filename", &outfile );
00043 
00044     opts.parseCommandLine( argc, argv );
00045 
00046     ErrorCode rval;
00047     Core* mb = new Core();
00048 
00049     rval = mb->load_file( inputfile.c_str() );MB_CHK_SET_ERR( rval, "can't load input file" );
00050 
00051     std::cout << " opened " << inputfile << " with initial h5m data.\n";
00052     // open the netcdf file, and see if it has that variable we are looking for
00053 
00054     Range nodes;
00055     rval = mb->get_entities_by_dimension( 0, 0, nodes );MB_CHK_SET_ERR( rval, "can't get nodes" );
00056 
00057     Range edges;
00058     rval = mb->get_entities_by_dimension( 0, 1, edges );MB_CHK_SET_ERR( rval, "can't get edges" );
00059 
00060     Range cells;
00061     rval = mb->get_entities_by_dimension( 0, 2, cells );MB_CHK_SET_ERR( rval, "can't get cells" );
00062 
00063     std::cout << " it has " << nodes.size() << " vertices " << edges.size() << " edges " << cells.size() << " cells\n";
00064 
00065     // construct maps between global id and handles
00066     std::map< int, EntityHandle > vGidHandle;
00067     std::map< int, EntityHandle > eGidHandle;
00068     std::map< int, EntityHandle > cGidHandle;
00069     std::vector< int > gids;
00070     Tag gid;
00071     rval = mb->tag_get_handle( "GLOBAL_ID", gid );MB_CHK_SET_ERR( rval, "can't get global id tag" );
00072     gids.resize( nodes.size() );
00073     rval = mb->tag_get_data( gid, nodes, &gids[0] );MB_CHK_SET_ERR( rval, "can't get global id on vertices" );
00074     int i = 0;
00075     for( Range::iterator vit = nodes.begin(); vit != nodes.end(); vit++ )
00076     {
00077         vGidHandle[gids[i++]] = *vit;
00078     }
00079 
00080     gids.resize( edges.size() );
00081     rval = mb->tag_get_data( gid, edges, &gids[0] );MB_CHK_SET_ERR( rval, "can't get global id on edges" );
00082     i = 0;
00083     for( Range::iterator vit = edges.begin(); vit != edges.end(); vit++ )
00084     {
00085         eGidHandle[gids[i++]] = *vit;
00086     }
00087 
00088     gids.resize( cells.size() );
00089     rval = mb->tag_get_data( gid, cells, &gids[0] );MB_CHK_SET_ERR( rval, "can't get global id on cells" );
00090     i = 0;
00091     for( Range::iterator vit = cells.begin(); vit != cells.end(); vit++ )
00092     {
00093         cGidHandle[gids[i++]] = *vit;
00094     }
00095 
00096     if( chunkfile_name.length() > 0 )
00097     {
00098 
00099         // Open chunk file
00100         ifstream inFile;
00101 
00102         inFile.open( chunkfile_name.c_str() );
00103         if( !inFile )
00104         {
00105             cout << "Unable to open chunk file";
00106             exit( 1 );  // terminate with error
00107         }
00108         Tag pTag, cTag;
00109         int def_val = -1;
00110         rval        = mb->tag_get_handle( "ProcID", 1, MB_TYPE_INTEGER, pTag, MB_TAG_CREAT | MB_TAG_DENSE, &def_val );MB_CHK_SET_ERR( rval, "can't define processor tag" );
00111         rval = mb->tag_get_handle( "ChunkID", 1, MB_TYPE_INTEGER, cTag, MB_TAG_CREAT | MB_TAG_DENSE, &def_val );MB_CHK_SET_ERR( rval, "can't define chunk tag" );
00112 
00113         int proc, lcid, ncols;
00114         while( inFile >> proc )
00115         {
00116             inFile >> lcid >> ncols;
00117             int Gid;
00118             for( i = 0; i < ncols; i++ )
00119             {
00120                 inFile >> Gid;
00121                 EntityHandle cell = cGidHandle[Gid];
00122                 rval              = mb->tag_set_data( pTag, &cell, 1, &proc );MB_CHK_SET_ERR( rval, "can't set proc tag" );
00123                 rval = mb->tag_set_data( cTag, &cell, 1, &lcid );MB_CHK_SET_ERR( rval, "can't set chunk tag" );
00124             }
00125         }
00126 
00127         inFile.close();
00128     }
00129 
00130     if( gsmapfile.length() > 0 )
00131     {
00132 
00133         // Open chunk file
00134         ifstream inFile;
00135 
00136         inFile.open( gsmapfile.c_str() );
00137         if( !inFile )
00138         {
00139             cout << "Unable to open gsmap file";
00140             exit( 1 );  // terminate with error
00141         }
00142         Tag pTag, cTag;
00143         int def_val             = -1;
00144         std::string procTagName = gsmapfile + "_proc";
00145         rval =
00146             mb->tag_get_handle( procTagName.c_str(), 1, MB_TYPE_INTEGER, pTag, MB_TAG_CREAT | MB_TAG_DENSE, &def_val );MB_CHK_SET_ERR( rval, "can't define processor tag" );
00147         std::string segTagName = gsmapfile + "_seg";
00148         rval =
00149             mb->tag_get_handle( segTagName.c_str(), 1, MB_TYPE_INTEGER, cTag, MB_TAG_CREAT | MB_TAG_DENSE, &def_val );MB_CHK_SET_ERR( rval, "can't define segment tag" );
00150 
00151         int compid, ngseg, gsize;
00152         inFile >> compid >> ngseg >> gsize;
00153         for( i = 1; i <= ngseg; i++ )
00154         {
00155             int start, len, pe;
00156             inFile >> start >> len >> pe;
00157             int Gid;
00158             for( int j = 0; j < len; j++ )
00159             {
00160                 Gid               = start + j;
00161                 EntityHandle cell = cGidHandle[Gid];
00162                 rval              = mb->tag_set_data( pTag, &cell, 1, &pe );MB_CHK_SET_ERR( rval, "can't set proc tag" );
00163                 rval = mb->tag_set_data( cTag, &cell, 1, &i );MB_CHK_SET_ERR( rval, "can't set segment tag" );
00164             }
00165         }
00166 
00167         inFile.close();
00168     }
00169 
00170     rval = mb->write_file( outfile.c_str() );MB_CHK_SET_ERR( rval, "can't write file" );
00171     std::cout << " wrote file " << outfile << "\n";
00172     return 0;
00173 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines