#include "IO/netcdf.h" #include "common/Utilities.h" #include "common/MPI_Helpers.h" #include "ProfilerApp.h" #ifdef USE_NETCDF #include #include #define CHECK_NC_ERR( ERR ) \ do { \ if ( ERR != NC_NOERR ) { \ std::string msg = "Error calling netcdf routine: "; \ msg += nc_strerror( ERR ); \ ERROR( msg ); \ } \ } while (0) namespace netcdf { // Convert nc_type to VariableType static inline VariableType convertType( nc_type type ) { VariableType type2 = UNKNOWN; if ( type == NC_BYTE ) type2 = BYTE; else if ( type == NC_CHAR ) type2 = STRING; else if ( type == NC_SHORT ) type2 = SHORT; else if ( type == NC_USHORT ) type2 = USHORT; else if ( type == NC_INT ) type2 = INT; else if ( type == NC_UINT ) type2 = UINT; else if ( type == NC_INT64 ) type2 = INT64; else if ( type == NC_UINT64 ) type2 = UINT64; else if ( type == NC_FLOAT ) type2 = FLOAT; else if ( type == NC_DOUBLE ) type2 = DOUBLE; else ERROR("Unknown type"); return type2; } // Get nc_type from the template template inline nc_type getType(); template<> inline nc_type getType() { return NC_CHAR; } template<> inline nc_type getType() { return NC_SHORT; } template<> inline nc_type getType() { return NC_INT; } template<> inline nc_type getType() { return NC_FLOAT; } template<> inline nc_type getType() { return NC_DOUBLE; } // Function to reverse an array template inline std::vector reverse( const std::vector& x ) { std::vector y(x.size()); for (size_t i=0; i inline std::vector convert( const std::vector& x ) { std::vector y(x.size()); for (size_t i=0; i(x[i]); return y; } /**************************************************** * Convert the VariableType to a string * ****************************************************/ std::string VariableTypeName( VariableType type ) { if ( type == BYTE ) return "BYTE"; else if ( type == SHORT ) return "SHORT"; else if ( type == USHORT ) return "USHORT"; else if ( type == INT ) return "INT"; else if ( type == UINT ) return "UINT"; else if ( type == INT64 ) return "INT64"; else if ( type == UINT64 ) return "UINT64"; else if ( type == FLOAT ) return "FLOAT"; else if ( type == DOUBLE ) return "DOUBLE"; else if ( type == STRING ) return "STRING"; return "Unknown"; } /**************************************************** * Open/close a file * ****************************************************/ int open( const std::string& filename, FileMode mode, MPI_Comm comm ) { int fid = 0; if ( comm == MPI_COMM_NULL ) { if ( mode == READ ) { int err = nc_open( filename.c_str(), NC_NOWRITE, &fid ); CHECK_NC_ERR( err ); } else if ( mode == WRITE ) { int err = nc_open( filename.c_str(), NC_WRITE, &fid ); CHECK_NC_ERR( err ); } else if ( mode == CREATE ) { int err = nc_create( filename.c_str(), NC_SHARE|NC_64BIT_OFFSET, &fid ); CHECK_NC_ERR( err ); } else { ERROR("Unknown file mode"); } } else { if ( mode == READ ) { int err = nc_open_par( filename.c_str(), NC_MPIPOSIX, comm, MPI_INFO_NULL, &fid ); CHECK_NC_ERR( err ); } else if ( mode == WRITE ) { int err = nc_open_par( filename.c_str(), NC_WRITE|NC_MPIPOSIX, comm, MPI_INFO_NULL, &fid ); CHECK_NC_ERR( err ); } else if ( mode == CREATE ) { int err = nc_create_par( filename.c_str(), NC_NETCDF4|NC_MPIIO, comm, MPI_INFO_NULL, &fid ); CHECK_NC_ERR( err ); } else { ERROR("Unknown file mode"); } } return fid; } void close( int fid ) { int err = nc_close( fid ); if ( err != NC_NOERR ) ERROR("Error closing file"); } /**************************************************** * Query basic properties * ****************************************************/ static std::vector getDimVar( int fid, int varid ) { int ndim = 0; int err = nc_inq_varndims( fid, varid, &ndim ); CHECK_NC_ERR( err ); std::vector dims(ndim,0); int dimid[64] = {-1}; err = nc_inq_vardimid( fid, varid, dimid ); CHECK_NC_ERR( err ); for (int i=0; i getVarDim( int fid, const std::string& var ) { return getDimVar( fid, getVarID( fid, var ) ); } std::vector getAttDim( int fid, const std::string& att ) { std::vector dim(1,0); int err = nc_inq_attlen( fid, NC_GLOBAL, att.c_str(), dim.data() ); CHECK_NC_ERR( err ); return dim; } std::vector getVarNames( int fid ) { int nvar; int err = nc_inq( fid, NULL, &nvar, NULL, NULL ); CHECK_NC_ERR( err ); std::vector vars(nvar); for (int i=0; i getAttNames( int fid ) { int natt; int err = nc_inq( fid, NULL, NULL, &natt, NULL ); CHECK_NC_ERR( err ); std::vector att(natt); for (int i=0; i Array getVar( int fid, const std::string& var ) { PROFILE_START("getVar"); Array x( reverse(getVarDim(fid,var)) ); int err = nc_get_var_ushort( fid, getVarID(fid,var), x.data() ); CHECK_NC_ERR( err ); PROFILE_STOP("getVar"); return x.reverseDim(); } template<> Array getVar( int fid, const std::string& var ) { PROFILE_START("getVar"); Array x( reverse(getVarDim(fid,var)) ); int err = nc_get_var_short( fid, getVarID(fid,var), x.data() ); CHECK_NC_ERR( err ); PROFILE_STOP("getVar"); return x.reverseDim(); } template<> Array getVar( int fid, const std::string& var ) { PROFILE_START("getVar"); Array x( reverse(getVarDim(fid,var)) ); int err = nc_get_var_uint( fid, getVarID(fid,var), x.data() ); CHECK_NC_ERR( err ); PROFILE_STOP("getVar"); return x.reverseDim(); } template<> Array getVar( int fid, const std::string& var ) { PROFILE_START("getVar"); Array x( reverse(getVarDim(fid,var)) ); int err = nc_get_var_int( fid, getVarID(fid,var), x.data() ); CHECK_NC_ERR( err ); PROFILE_STOP("getVar"); return x.reverseDim(); } template<> Array getVar( int fid, const std::string& var ) { PROFILE_START("getVar"); Array x( reverse(getVarDim(fid,var)) ); int err = nc_get_var_float( fid, getVarID(fid,var), x.data() ); CHECK_NC_ERR( err ); PROFILE_STOP("getVar"); return x.reverseDim(); } template<> Array getVar( int fid, const std::string& var ) { PROFILE_START("getVar"); Array x( reverse(getVarDim(fid,var)) ); int err = nc_get_var_double( fid, getVarID(fid,var), x.data() ); CHECK_NC_ERR( err ); PROFILE_STOP("getVar"); return x.reverseDim(); } template<> Array getVar( int fid, const std::string& var ) { PROFILE_START("getVar"); Array x( reverse(getVarDim(fid,var)) ); int err = nc_get_var_text( fid, getVarID(fid,var), x.data() ); CHECK_NC_ERR( err ); PROFILE_STOP("getVar"); return x.reverseDim(); } template<> Array getVar( int fid, const std::string& var ) { PROFILE_START("getVar"); Array tmp = getVar( fid, var ); std::vector dim = {tmp.size(0), tmp.size(1), tmp.size(2) }; if ( dim.size() == 1 ) dim[0] = 1; else dim.erase( dim.begin() ); Array text(dim); for (size_t i=0; i"); return text; } static inline void get_stride_args( const std::vector& start, const std::vector& count, const std::vector& stride, size_t *startp, size_t *countp, ptrdiff_t *stridep ) { for (size_t i=0; i int nc_get_vars_TYPE( int fid, int varid, const size_t start[], const size_t count[], const ptrdiff_t stride[], TYPE *ptr ); template<> int nc_get_vars_TYPE( int fid, int varid, const size_t start[], const size_t count[], const ptrdiff_t stride[], short *ptr ) { return nc_get_vars_short( fid, varid, start, count, stride, ptr ); } template<> int nc_get_vars_TYPE( int fid, int varid, const size_t start[], const size_t count[], const ptrdiff_t stride[], int *ptr ) { return nc_get_vars_int( fid, varid, start, count, stride, ptr ); } template<> int nc_get_vars_TYPE( int fid, int varid, const size_t start[], const size_t count[], const ptrdiff_t stride[], float *ptr ) { return nc_get_vars_float( fid, varid, start, count, stride, ptr ); } template<> int nc_get_vars_TYPE( int fid, int varid, const size_t start[], const size_t count[], const ptrdiff_t stride[], double *ptr ) { return nc_get_vars_double( fid, varid, start, count, stride, ptr ); } template Array getVar( int fid, const std::string& var, const std::vector& start, const std::vector& count, const std::vector& stride ) { PROFILE_START("getVar<> (strided)"); std::vector var_size = getVarDim( fid, var ); for (int d=0; d<(int)var_size.size(); d++) { if ( start[d]<0 || start[d]+stride[d]*(count[d]-1)>(int)var_size[d] ) { int rank = comm_rank(MPI_COMM_WORLD); char tmp[1000]; sprintf(tmp,"%i: Range exceeded array dimension:\n" " start[%i]=%i, count[%i]=%i, stride[%i]=%i, var_size[%i]=%i", rank,d,start[d],d,count[d],d,stride[d],d,(int)var_size[d]); ERROR(tmp); } } Array x( reverse(convert(count)) ); size_t startp[10], countp[10]; ptrdiff_t stridep[10]; get_stride_args( start, count, stride, startp, countp, stridep ); int err = nc_get_vars_TYPE( fid, getVarID(fid,var), startp, countp, stridep, x.data() ); CHECK_NC_ERR( err ); PROFILE_STOP("getVar<> (strided)"); return x.reverseDim(); } template Array getVar( int, const std::string&, const std::vector&, const std::vector&, const std::vector& ); template Array getVar( int, const std::string&, const std::vector&, const std::vector&, const std::vector& ); template Array getVar( int, const std::string&, const std::vector&, const std::vector&, const std::vector& ); template Array getVar( int, const std::string&, const std::vector&, const std::vector&, const std::vector& ); /**************************************************** * Read an attribute * ****************************************************/ template<> Array getAtt( int fid, const std::string& att ) { PROFILE_START("getAtt"); Array x( getAttDim(fid,att) ); int err = nc_get_att_double( fid, NC_GLOBAL, att.c_str(), x.data() ); CHECK_NC_ERR( err ); PROFILE_STOP("getAtt"); return x; } template<> Array getAtt( int fid, const std::string& att ) { PROFILE_START("getAtt"); char *tmp = new char[getAttDim(fid,att)[0]]; Array x(1); x(0) = tmp; delete [] tmp; PROFILE_STOP("getAtt"); return x; } /**************************************************** * Write an array to a file * ****************************************************/ std::vector defDim( int fid, const std::vector& names, const std::vector& dims ) { std::vector dimid(names.size(),0); for (size_t i=0; i void write( int fid, const std::string& var, const std::vector& dimids, const Array& data, const RankInfoStruct& info ) { // Define the variable int varid = 0; int err = nc_def_var( fid, var.c_str(), getType(), data.ndim(), dimids.data(), &varid ); CHECK_NC_ERR( err ); // exit define mode err = nc_enddef( fid ); CHECK_NC_ERR( err ); // set the access method to use MPI/PnetCDF collective I/O err = nc_var_par_access( fid, varid, NC_INDEPENDENT ); CHECK_NC_ERR( err ); // parallel write: each process writes its subarray to the file auto x = data.reverseDim(); std::vector count = { data.size(0), data.size(1), data.size(2) }; std::vector start = { info.ix*data.size(0), info.jy*data.size(1), info.kz*data.size(2) }; nc_put_vara( fid, varid, start.data(), count.data(), x.data() ); } template void write( int fid, const std::string& var, const std::vector& dimids, const Array& data, const RankInfoStruct& info ); template void write( int fid, const std::string& var, const std::vector& dimids, const Array& data, const RankInfoStruct& info ); template void write( int fid, const std::string& var, const std::vector& dimids, const Array& data, const RankInfoStruct& info ); template void write( int fid, const std::string& var, const std::vector& dimids, const Array& data, const RankInfoStruct& info ); }; // netcdf namespace #else #endif