MOAB: Mesh Oriented datABase
(version 5.4.1)
|
#include "moab/MOABConfig.h"
#include "moab/ProgOptions.hpp"
#include "netcdf.h"
#include <cmath>
#include <iomanip>
#include <Eigen/Sparse>
Go to the source code of this file.
Defines | |
#define | ERR_NC(e) |
#define | INS_ID(stringvar, prefix, id) sprintf( stringvar, prefix, id ) |
#define | GET_DIM1(ncdim, name, val) |
#define | GET_DIMB1(ncdim, name, varname, id, val) |
#define | GET_VAR1(name, id, dims) |
#define | GET_1D_INT_VAR1(name, id, vals) |
#define | GET_1D_DBL_VAR1(name, id, vals) |
#define | GET_DIM2(ncdim, name, val) |
#define | GET_DIMB2(ncdim, name, varname, id, val) |
#define | GET_VAR2(name, id, dims) |
#define | GET_1D_INT_VAR2(name, id, vals) |
#define | GET_1D_DBL_VAR2(name, id, vals) |
#define | GET_2D_DBL_VAR1(name, id, vals) |
#define | GET_2D_DBL_VAR2(name, id, vals) |
Typedefs | |
typedef Eigen::Map < Eigen::VectorXd > | EigenV |
Functions | |
void | diff_vect (const char *var_name, int n) |
void | diff_2d_vect (const char *var_name, int n) |
int | main (int argc, char *argv[]) |
Variables | |
int | ncFile1 |
int | ncFile2 |
#define ERR_NC | ( | e | ) |
{ \
printf( "Error: %s\n", nc_strerror( e ) ); \
exit( 2 ); \
}
Definition at line 27 of file compareMaps.cpp.
Referenced by main().
#define GET_1D_DBL_VAR1 | ( | name, | |
id, | |||
vals | |||
) |
{ \ std::vector< int > dum_dims; \ GET_VAR1( name, id, dum_dims ); \ if( -1 != ( id ) ) \ { \ size_t ntmp; \ int dvfail = nc_inq_dimlen( ncFile1, dum_dims[0], &ntmp ); \ if( NC_NOERR != dvfail ) \ { \ ERR_NC( dvfail ) \ } \ ( vals ).resize( ntmp ); \ size_t ntmp1 = 0; \ dvfail = nc_get_vara_double( ncFile1, id, &ntmp1, &ntmp, &( vals )[0] ); \ if( NC_NOERR != dvfail ) \ { \ ERR_NC( dvfail ) \ } \ } \ }
Definition at line 103 of file compareMaps.cpp.
Referenced by diff_vect(), and main().
#define GET_1D_DBL_VAR2 | ( | name, | |
id, | |||
vals | |||
) |
{ \ std::vector< int > dum_dims; \ GET_VAR2( name, id, dum_dims ); \ if( -1 != ( id ) ) \ { \ size_t ntmp; \ int dvfail = nc_inq_dimlen( ncFile2, dum_dims[0], &ntmp ); \ if( NC_NOERR != dvfail ) \ { \ ERR_NC( dvfail ) \ } \ ( vals ).resize( ntmp ); \ size_t ntmp1 = 0; \ dvfail = nc_get_vara_double( ncFile2, id, &ntmp1, &ntmp, &( vals )[0] ); \ if( NC_NOERR != dvfail ) \ { \ ERR_NC( dvfail ) \ } \ } \ }
Definition at line 190 of file compareMaps.cpp.
Referenced by diff_vect(), and main().
#define GET_1D_INT_VAR1 | ( | name, | |
id, | |||
vals | |||
) |
{ \ GET_VAR1( name, id, vals ); \ if( -1 != ( id ) ) \ { \ size_t ntmp; \ int ivfail = nc_inq_dimlen( ncFile1, ( vals )[0], &ntmp ); \ if( NC_NOERR != ivfail ) \ { \ ERR_NC( ivfail ) \ } \ ( vals ).resize( ntmp ); \ size_t ntmp1 = 0; \ ivfail = nc_get_vara_int( ncFile1, id, &ntmp1, &ntmp, &( vals )[0] ); \ if( NC_NOERR != ivfail ) \ { \ ERR_NC( ivfail ) \ } \ } \ }
Definition at line 82 of file compareMaps.cpp.
Referenced by main().
#define GET_1D_INT_VAR2 | ( | name, | |
id, | |||
vals | |||
) |
{ \ GET_VAR2( name, id, vals ); \ if( -1 != ( id ) ) \ { \ size_t ntmp; \ int ivfail = nc_inq_dimlen( ncFile2, ( vals )[0], &ntmp ); \ if( NC_NOERR != ivfail ) \ { \ ERR_NC( ivfail ) \ } \ ( vals ).resize( ntmp ); \ size_t ntmp1 = 0; \ ivfail = nc_get_vara_int( ncFile2, id, &ntmp1, &ntmp, &( vals )[0] ); \ if( NC_NOERR != ivfail ) \ { \ ERR_NC( ivfail ) \ } \ } \ }
Definition at line 169 of file compareMaps.cpp.
Referenced by main().
#define GET_2D_DBL_VAR1 | ( | name, | |
id, | |||
vals | |||
) |
{ \ std::vector< int > dum_dims; \ GET_VAR1( name, id, dum_dims ); \ if( -1 != ( id ) ) \ { \ size_t ntmp[2]; \ int dvfail = nc_inq_dimlen( ncFile1, dum_dims[0], &ntmp[0] ); \ if( NC_NOERR != dvfail ) \ { \ ERR_NC( dvfail ) \ } \ dvfail = nc_inq_dimlen( ncFile1, dum_dims[1], &ntmp[1] ); \ if( NC_NOERR != dvfail ) \ { \ ERR_NC( dvfail ) \ } \ ( vals ).resize( ntmp[0] * ntmp[1] ); \ size_t ntmp1[2] = { 0, 0 }; \ dvfail = nc_get_vara_double( ncFile1, id, ntmp1, ntmp, &( vals )[0] ); \ if( NC_NOERR != dvfail ) \ { \ ERR_NC( dvfail ) \ } \ } \ }
Definition at line 212 of file compareMaps.cpp.
Referenced by diff_2d_vect().
#define GET_2D_DBL_VAR2 | ( | name, | |
id, | |||
vals | |||
) |
{ \ std::vector< int > dum_dims; \ GET_VAR2( name, id, dum_dims ); \ if( -1 != ( id ) ) \ { \ size_t ntmp[2]; \ int dvfail = nc_inq_dimlen( ncFile2, dum_dims[0], &ntmp[0] ); \ if( NC_NOERR != dvfail ) \ { \ ERR_NC( dvfail ) \ } \ dvfail = nc_inq_dimlen( ncFile2, dum_dims[1], &ntmp[1] ); \ if( NC_NOERR != dvfail ) \ { \ ERR_NC( dvfail ) \ } \ ( vals ).resize( ntmp[0] * ntmp[1] ); \ size_t ntmp1[2] = { 0, 0 }; \ dvfail = nc_get_vara_double( ncFile2, id, ntmp1, ntmp, &( vals )[0] ); \ if( NC_NOERR != dvfail ) \ { \ ERR_NC( dvfail ) \ } \ } \ }
Definition at line 239 of file compareMaps.cpp.
Referenced by diff_2d_vect().
#define GET_DIM1 | ( | ncdim, | |
name, | |||
val | |||
) |
{ \ int gdfail = nc_inq_dimid( ncFile1, name, &( ncdim ) ); \ if( NC_NOERR == gdfail ) \ { \ size_t tmp_val; \ gdfail = nc_inq_dimlen( ncFile1, ncdim, &tmp_val ); \ if( NC_NOERR != gdfail ) \ { \ ERR_NC( gdfail ) \ } \ else \ ( val ) = tmp_val; \ } \ else \ ( val ) = 0; \ }
Definition at line 40 of file compareMaps.cpp.
Referenced by main().
#define GET_DIM2 | ( | ncdim, | |
name, | |||
val | |||
) |
{ \ int gdfail = nc_inq_dimid( ncFile2, name, &( ncdim ) ); \ if( NC_NOERR == gdfail ) \ { \ size_t tmp_val; \ gdfail = nc_inq_dimlen( ncFile2, ncdim, &tmp_val ); \ if( NC_NOERR != gdfail ) \ { \ ERR_NC( gdfail ) \ } \ else \ ( val ) = tmp_val; \ } \ else \ ( val ) = 0; \ }
Definition at line 127 of file compareMaps.cpp.
Referenced by main().
#define GET_DIMB1 | ( | ncdim, | |
name, | |||
varname, | |||
id, | |||
val | |||
) |
INS_ID( name, varname, id ); \ GET_DIM1( ncdim, name, val );
Definition at line 58 of file compareMaps.cpp.
#define GET_DIMB2 | ( | ncdim, | |
name, | |||
varname, | |||
id, | |||
val | |||
) |
INS_ID( name, varname, id ); \ GET_DIM2( ncdim, name, val );
Definition at line 145 of file compareMaps.cpp.
#define GET_VAR1 | ( | name, | |
id, | |||
dims | |||
) |
{ \ ( id ) = -1; \ int gvfail = nc_inq_varid( ncFile1, name, &( id ) ); \ if( NC_NOERR == gvfail ) \ { \ int ndims; \ gvfail = nc_inq_varndims( ncFile1, id, &ndims ); \ if( NC_NOERR == gvfail ) \ { \ ( dims ).resize( ndims ); \ gvfail = nc_inq_vardimid( ncFile1, id, &( dims )[0] ); \ if( NC_NOERR != gvfail ) \ { \ ERR_NC( gvfail ) \ } \ } \ } \ }
Definition at line 62 of file compareMaps.cpp.
#define GET_VAR2 | ( | name, | |
id, | |||
dims | |||
) |
{ \ ( id ) = -1; \ int gvfail = nc_inq_varid( ncFile2, name, &( id ) ); \ if( NC_NOERR == gvfail ) \ { \ int ndims; \ gvfail = nc_inq_varndims( ncFile2, id, &ndims ); \ if( NC_NOERR == gvfail ) \ { \ ( dims ).resize( ndims ); \ gvfail = nc_inq_vardimid( ncFile2, id, &( dims )[0] ); \ if( NC_NOERR != gvfail ) \ { \ ERR_NC( gvfail ) \ } \ } \ } \ }
Definition at line 149 of file compareMaps.cpp.
#define INS_ID | ( | stringvar, | |
prefix, | |||
id | |||
) | sprintf( stringvar, prefix, id ) |
Definition at line 38 of file compareMaps.cpp.
typedef Eigen::Map< Eigen::VectorXd > EigenV |
Definition at line 266 of file compareMaps.cpp.
void diff_2d_vect | ( | const char * | var_name, |
int | n | ||
) |
Definition at line 289 of file compareMaps.cpp.
References GET_2D_DBL_VAR1, and GET_2D_DBL_VAR2.
Referenced by main().
{ // compare frac_a between maps std::vector< double > fraca1( n ), fraca2( n ); int idfa1, idfa2; GET_2D_DBL_VAR1( var_name, idfa1, fraca1 ); EigenV fa1( fraca1.data(), n ); GET_2D_DBL_VAR2( var_name, idfa2, fraca2 ); EigenV fa2( fraca2.data(), n ); std::cout << var_name << " diff norm: " << ( fa1 - fa2 ).norm() << "\n"; return; }
void diff_vect | ( | const char * | var_name, |
int | n | ||
) |
Definition at line 268 of file compareMaps.cpp.
References GET_1D_DBL_VAR1, and GET_1D_DBL_VAR2.
Referenced by main().
{ // compare frac_a between maps std::vector< double > fraca1( n ), fraca2( n ); int idfa1, idfa2; GET_1D_DBL_VAR1( var_name, idfa1, fraca1 ); EigenV fa1( fraca1.data(), n ); GET_1D_DBL_VAR2( var_name, idfa2, fraca2 ); EigenV fa2( fraca2.data(), n ); EigenV diff( fraca2.data(), n ); diff = fa1 - fa2; int imin, imax; double minV = diff.minCoeff( &imin ); double maxV = diff.maxCoeff( &imax ); std::cout << var_name << " diff norm: " << diff.norm() << " min at " << imin << " : " << minV << " max at " << imax << " : " << maxV << "\n"; return; }
int main | ( | int | argc, |
char * | argv[] | ||
) |
Definition at line 301 of file compareMaps.cpp.
References ProgOptions::addOpt(), diff_2d_vect(), diff_vect(), ERR_NC, fail(), GET_1D_DBL_VAR1, GET_1D_DBL_VAR2, GET_1D_INT_VAR1, GET_1D_INT_VAR2, GET_DIM1, GET_DIM2, ncFile1, ncFile2, and ProgOptions::parseCommandLine().
{ ProgOptions opts; std::string inputfile1, inputfile2; opts.addOpt< std::string >( "firstMap,i", "input filename 1", &inputfile1 ); opts.addOpt< std::string >( "secondMap,j", "input second map", &inputfile2 ); int print_diff = 20; opts.addOpt< int >( "print_differences,p", "print differences ", &print_diff ); double fraction = 0.9; opts.addOpt< double >( "fraction_diff,f", "fraction threshold", &fraction ); opts.parseCommandLine( argc, argv ); // Open netcdf/exodus file int fail = nc_open( inputfile1.c_str(), 0, &ncFile1 ); if( NC_NOWRITE != fail ) { ERR_NC( fail ) } std::cout << " opened " << inputfile1 << " for map 1 \n"; fail = nc_open( inputfile2.c_str(), 0, &ncFile2 ); if( NC_NOWRITE != fail ) { ERR_NC( fail ) } std::cout << " opened " << inputfile2 << " for map 2 \n"; int temp_dim; int na1, nb1, ns1, na2, nb2, ns2, nv_a, nv_b, nv_a2, nv_b2; GET_DIM1( temp_dim, "n_a", na1 ); GET_DIM2( temp_dim, "n_a", na2 ); GET_DIM1( temp_dim, "n_b", nb1 ); GET_DIM2( temp_dim, "n_b", nb2 ); GET_DIM1( temp_dim, "n_s", ns1 ); GET_DIM2( temp_dim, "n_s", ns2 ); GET_DIM1( temp_dim, "nv_a", nv_a ); GET_DIM2( temp_dim, "nv_a", nv_a2 ); GET_DIM1( temp_dim, "nv_b", nv_b ); GET_DIM2( temp_dim, "nv_b", nv_b2 ); if( nv_a != nv_a2 || nv_b != nv_b2 ) { std::cout << " different nv dimensions:" << nv_a << " == " << nv_a2 << " or " << nv_b << " == " << nv_b2 << " bail out \n"; return 1; } std::cout << " n_a, n_b, n_s : " << na1 << ", " << nb1 << ", " << ns1 << " for map 1 \n"; std::cout << " n_a, n_b, n_s : " << na2 << ", " << nb2 << ", " << ns2 << " for map 2 \n"; if( na1 != na2 || nb1 != nb2 ) { std::cout << " different dimensions bail out \n"; return 1; } std::vector< int > col1( ns1 ), row1( ns1 ); std::vector< int > col2( ns2 ), row2( ns2 ); std::vector< double > val1( ns1 ), val2( ns2 ); int idrow1, idcol1, idrow2, idcol2, ids1, ids2; GET_1D_INT_VAR1( "row", idrow1, row1 ); GET_1D_INT_VAR1( "col", idcol1, col1 ); GET_1D_DBL_VAR1( "S", ids1, val1 ); GET_1D_INT_VAR2( "row", idrow2, row2 ); GET_1D_INT_VAR2( "col", idcol2, col2 ); GET_1D_DBL_VAR2( "S", ids2, val2 ); // first matrix typedef Eigen::Triplet< double > Triplet; std::vector< Triplet > tripletList; tripletList.reserve( ns1 ); for( int iv = 0; iv < ns1; iv++ ) { tripletList.push_back( Triplet( row1[iv] - 1, col1[iv] - 1, val1[iv] ) ); } Eigen::SparseMatrix< double > weight1( nb1, na1 ); weight1.setFromTriplets( tripletList.begin(), tripletList.end() ); weight1.makeCompressed(); if( ns1 != ns2 ) tripletList.resize( ns2 ); for( int iv = 0; iv < ns2; iv++ ) { tripletList[iv] = Triplet( row2[iv] - 1, col2[iv] - 1, val2[iv] ); } Eigen::SparseMatrix< double > weight2( nb2, na2 ); weight2.setFromTriplets( tripletList.begin(), tripletList.end() ); weight2.makeCompressed(); // default storage type is column major Eigen::SparseMatrix< double > diff = weight1 - weight2; diff.makeCompressed(); // is it needed or not ? auto coeffs = diff.coeffs(); double maxv = coeffs.maxCoeff(); double minv = coeffs.minCoeff(); std::cout << " euclidian norm for difference: " << diff.norm() << " \n squared norm for difference: " << diff.squaredNorm() << "\n" << " minv: " << minv << " maxv: " << maxv << "\n"; // print out the first 10 positions for which the value is outside 90% of min/max values double min_threshold = fraction * minv; double max_threshold = fraction * maxv; int counter = 0; std::cout << std::setprecision( 9 ); for( int k = 0; ( k < diff.outerSize() ); ++k ) // this is by column { for( Eigen::SparseMatrix< double >::InnerIterator it( diff, k ); ( it ) && ( counter < print_diff ); ++it ) { double val = it.value(); if( val <= min_threshold || val >= max_threshold ) { int row = it.row(); int col = it.col(); std::cout << " counter:" << counter << "\t col: " << col + 1 << "\t row: " << row + 1 << "\t diff: " << val; std::cout << "\t map1: " << weight1.coeffRef( row, col ) << "\t map2: " << weight2.coeffRef( row, col ) << "\n"; // row index counter++; } } } // compare frac_a between maps diff_vect( "frac_a", na1 ); diff_vect( "frac_b", nb1 ); diff_vect( "area_a", na1 ); diff_vect( "area_b", nb1 ); diff_vect( "yc_a", na1 ); diff_vect( "yc_b", nb1 ); diff_vect( "xc_a", na1 ); diff_vect( "xc_b", nb1 ); diff_2d_vect( "xv_a", na1 * nv_a ); diff_2d_vect( "yv_a", na1 * nv_a ); diff_2d_vect( "xv_b", nb1 * nv_b ); diff_2d_vect( "yv_b", nb1 * nv_b ); return 0; }
int ncFile1 |
Definition at line 36 of file compareMaps.cpp.
Referenced by main().
int ncFile2 |
Definition at line 125 of file compareMaps.cpp.
Referenced by main().