[MOAB-dev] commit/MOAB: danwu: Updated NC writer to transpose tag values of non-set variables (with timesteps, e.g. variable T) before writing.
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Fri Apr 4 16:29:01 CDT 2014
1 new commit in MOAB:
https://bitbucket.org/fathomteam/moab/commits/2bdfc24fa3c8/
Changeset: 2bdfc24fa3c8
Branch: ncwriter
User: danwu
Date: 2014-04-04 23:28:46
Summary: Updated NC writer to transpose tag values of non-set variables (with timesteps, e.g. variable T) before writing.
Affected #: 2 files
diff --git a/src/io/WriteNC.cpp b/src/io/WriteNC.cpp
index 30601c0..b7c0450 100644
--- a/src/io/WriteNC.cpp
+++ b/src/io/WriteNC.cpp
@@ -277,7 +277,7 @@ ErrorCode WriteNC::process_conventional_tags(EntityHandle fileSet)
rval = mbImpl->tag_get_by_ptr(tagMeshType, &fileSet, 1, &data, &sz);
ERRORR(rval, "Trouble getting values for conventional tag __MESH_TYPE.");
- p = static_cast<const char *>(data);
+ p = static_cast<const char*>(data);
grid_type = std::string(&p[0], sz);
dbgOut.tprintf(2, "mesh type: %s \n", grid_type.c_str());
@@ -596,7 +596,7 @@ ErrorCode WriteNC::collect_variable_data(std::vector<std::string>& var_names, st
Tag coordtag = 0;
rval = mbImpl->tag_get_handle(coordName.c_str(), coordtag);
ERRORR(rval, "Can't find one tag.");
- varCoordData.varTags.push_back(coordtag); // really, only one for these
+ varCoordData.varTags.push_back(coordtag); // Really, only one for these
const void* data;
int sizeCoordinate;
@@ -684,12 +684,12 @@ ErrorCode WriteNC::initialize_file(std::vector<std::string>& var_names)
ERRORR(MB_FAILURE, "Can't find coordinate variable requested.");
VarData& coordData = vit2->second;
- variableData.varDims[j] = coordData.varDims[0]; // this one, being a coordinate, is the only one
+ variableData.varDims[j] = coordData.varDims[0]; // This one, being a coordinate, is the only one
dbgOut.tprintf(2, " dimension with index %d name %s has ID %d \n",
j, dimName.c_str(), variableData.varDims[j]);
- variableData.writeStarts.push_back(0); // assume we will write all, so start at 0 for all dimensions
- variableData.writeCounts.push_back(coordData.sz); // again, write all; times will be one at a time
+ variableData.writeStarts.push_back(0); // Assume we will write all, so start at 0 for all dimensions
+ variableData.writeCounts.push_back(coordData.sz); // Again, write all; times will be one at a time
}
// Define the variable now:
@@ -785,7 +785,8 @@ ErrorCode WriteNC::write_values(std::vector<std::string>& var_names, EntityHandl
// Now look at requested var_names; if they have time, we will have a list, and write one at a time
// We may also need to gather, and transpose stuff
Range ents2d;
- // Get all entities of dimension 2 from set; assume now location is on cells;
+ // Get all entities of dimension 2 from set
+ // FIXME: assume now location is on cells (for variable T of CAM-EUL it is the case)
// Need to reorder stuff in the order from the file, also transpose from lev dimension
ErrorCode rval = mbImpl->get_entities_by_dimension(fileSet, 2, ents2d);
ERRORR(rval, "Can't get entities for 2d.");
@@ -800,14 +801,18 @@ ErrorCode WriteNC::write_values(std::vector<std::string>& var_names, EntityHandl
VarData& variableData = vit->second;
int numTimeSteps = (int)variableData.varTags.size();
if (variableData.has_tsteps) {
- variableData.writeCounts[0] = 1; // we will write one time step
+ // FIXME: assume now the variable has 4 dimensions as (time, lev, lat, lon)
+ // At each timestep, we need to transpose tag format (lat, lon, lev) back
+ // to NC format (lev, lat, lon) for writing
+ size_t ni = variableData.writeCounts[3]; // lon
+ size_t nj = variableData.writeCounts[2]; // lat
+ size_t nk = variableData.writeCounts[1]; // lev
+
+ variableData.writeCounts[0] = 1; // We will write one time step
for (int j = 0; j < numTimeSteps; j++) {
// We will write one time step, and count will be one; start will be different
// We will write values directly from tag_iterate, but we should also transpose for level
- // So that means deep copy for transpose
- // !!!!!!!!!!!!!!
- //
- // FIXME !!!!!!!!!!!
+ // so that means deep copy for transpose
variableData.writeStarts[0] = j; // This is time, again
int count;
void* dataptr;
@@ -815,21 +820,29 @@ ErrorCode WriteNC::write_values(std::vector<std::string>& var_names, EntityHandl
assert(count == (int)ents2d.size());
// Now write from memory directly
- // FIXME: we need to transpose and gather for multiple processors
+ // FIXME: we need to gather for multiple processors
int success = 0;
switch (variableData.varDataType) {
- case NC_DOUBLE:
+ case NC_DOUBLE: {
+ std::vector<double> tmpdoubledata(ni*nj*nk);
+ // Transpose (lat, lon, lev) back to (lev, lat, lon)
+ jik_to_kji(ni, nj, nk, &tmpdoubledata[0], (double*)(dataptr));
success = NCFUNCAP(_vara_double)(fileId, variableData.varId,
&variableData.writeStarts[0], &variableData.writeCounts[0],
- (double*)(dataptr));
+ &tmpdoubledata[0]);
ERRORS(success, "Failed to write double data.");
break;
- case NC_INT:
+ }
+ case NC_INT: {
+ std::vector<int> tmpintdata(ni*nj*nk);
+ // Transpose (lat, lon, lev) back to (lev, lat, lon)
+ jik_to_kji(ni, nj, nk, &tmpintdata[0], (int*)(dataptr));
success = NCFUNCAP(_vara_int)(fileId, variableData.varId,
&variableData.writeStarts[0], &variableData.writeCounts[0],
- (int*)(dataptr));
+ &tmpintdata[0]);
ERRORS(success, "Failed to write int data.");
break;
+ }
default:
success = 1;
break;
diff --git a/src/io/WriteNC.hpp b/src/io/WriteNC.hpp
index b027bf6..36207d7 100644
--- a/src/io/WriteNC.hpp
+++ b/src/io/WriteNC.hpp
@@ -16,9 +16,8 @@
#ifndef WRITENC_HPP_
#define WRITENC_HPP_
-
#ifndef IS_BUILDING_MB
-//#error "ReadNC.hpp isn't supposed to be included into an application"
+//#error "WriteNC.hpp isn't supposed to be included into an application"
#endif
#include <vector>
@@ -39,26 +38,15 @@
#include "pnetcdf.h"
#define NCFUNC(func) ncmpi_ ## func
-//! Collective I/O mode get
-#define NCFUNCAG(func) ncmpi_get ## func ## _all
-
//! Collective I/O mode put
#define NCFUNCAP(func) ncmpi_put ## func ## _all
-//! Independent I/O mode get
-#define NCFUNCG(func) ncmpi_get ## func
-
-//! Nonblocking get (request aggregation), used so far only for ucd mesh
-#define NCFUNCREQG(func) ncmpi_iget ## func
-
#define NCDF_SIZE MPI_Offset
#define NCDF_DIFF MPI_Offset
#else
#include "netcdf.h"
#define NCFUNC(func) nc_ ## func
-#define NCFUNCAG(func) nc_get ## func
#define NCFUNCAP(func) nc_put ## func
-#define NCFUNCG(func) nc_get ## func
#define NCDF_SIZE size_t
#define NCDF_DIFF ptrdiff_t
#endif
@@ -72,7 +60,6 @@ class NCWriteHelper;
/**
* \brief Export NC files.
*/
-
class WriteNC : public WriterIface
{
friend class NCWriteHelper;
@@ -80,36 +67,36 @@ class WriteNC : public WriterIface
public:
- //! factory method
- static WriterIface* factory( Interface* );
+ //! Factory method
+ static WriterIface* factory(Interface*);
- //! Constructor
- WriteNC(Interface *impl);
+ //! Constructor
+ WriteNC(Interface* impl = NULL);
- //! Destructor
+ //! Destructor
virtual ~WriteNC();
-
- //! writes out a file
- ErrorCode write_file(const char *file_name,
- const bool overwrite,
- const FileOptions& opts,
- const EntityHandle *output_list,
- const int num_sets,
- const std::vector<std::string>& qa_list,
- const Tag* tag_list = NULL,
- int num_tags = 0,
- int export_dimension = 3);
-
+ //! Writes out a file
+ ErrorCode write_file(const char* file_name,
+ const bool overwrite,
+ const FileOptions& opts,
+ const EntityHandle* output_list,
+ const int num_sets,
+ const std::vector<std::string>& qa_list,
+ const Tag* tag_list = NULL,
+ int num_tags = 0,
+ int export_dimension = 3);
private:
+ //! ENTLOCNSEDGE for north/south edge
+ //! ENTLOCWEEDGE for west/east edge
enum EntityLocation {ENTLOCVERT = 0, ENTLOCNSEDGE, ENTLOCEWEDGE, ENTLOCFACE, ENTLOCSET, ENTLOCEDGE, ENTLOCREGION};
class AttData
{
public:
- AttData() : attId(-1), attLen(0), attVarId(-2) {}
+ AttData() : attId(-1), attLen(0), attVarId(-2), attDataType(NC_NAT) {}
int attId;
NCDF_SIZE attLen;
int attVarId;
@@ -127,23 +114,23 @@ private:
std::map<std::string, AttData> varAtts;
std::string varName;
std::vector<Tag> varTags; // Tags created for this variable, e.g. one tag per timestep
- std::vector<void*> memoryHogs; // these will point to the real data; fill before writing the data
- std::vector<NCDF_SIZE> writeStarts; // Starting index for reading data values along each dimension
- std::vector<NCDF_SIZE> writeCounts; // Number of data values to be read along each dimension
+ std::vector<void*> memoryHogs; // These will point to the real data; fill before writing the data
+ std::vector<NCDF_SIZE> writeStarts; // Starting index for writing data values along each dimension
+ std::vector<NCDF_SIZE> writeCounts; // Number of data values to be written along each dimension
int entLoc;
int numLev;
int sz;
bool has_tsteps; // Indicate whether timestep numbers are appended to tag names
};
- // this info will be reconstructed from metadata stored on conventional fileSet tags
+ // This info will be reconstructed from metadata stored on conventional fileSet tags
//! Dimension names
std::vector<std::string> dimNames;
//! Dimension lengths
std::vector<int> dimLens;
- // will collect used dimensions (ccordinate variables)
+ // Will collect used dimensions (coordinate variables)
std::set<std::string> usedCoordinates;
//! Global attribs
@@ -152,39 +139,51 @@ private:
//! Variable info
std::map<std::string, VarData> varInfo;
- ErrorCode parse_options(const FileOptions& opts, std::vector<std::string>& var_names, std::vector<int>& tstep_nums,
- std::vector<double>& tstep_vals);
+ ErrorCode parse_options(const FileOptions& opts,
+ std::vector<std::string>& var_names,
+ std::vector<int>& tstep_nums,
+ std::vector<double>& tstep_vals);
/*
- * map out the header, from tags on file set; it is the inverse process from
+ * Map out the header, from tags on file set; it is the inverse process from
* ErrorCode NCHelper::create_conventional_tags
*/
ErrorCode process_conventional_tags(EntityHandle fileSet);
- ErrorCode process_concatenated_attribute(const void * gattptr, int globalAttSz, std::vector<int> & gattLen,
- std::map<std::string, AttData> & globalAtts);
+ ErrorCode process_concatenated_attribute(const void* gattptr, int globalAttSz,
+ std::vector<int>& gattLen,
+ std::map<std::string, AttData>& globalAtts);
- // will collect data; it should be only on gather processor, but for the time being, collect
+ // Will collect data; it should be only on gather processor, but for the time being, collect
// for everybody
- ErrorCode collect_variable_data( std::vector<std::string>& var_names, std::vector<int>& tstep_nums,
- std::vector<double>& tstep_vals, EntityHandle fileSet);
+ ErrorCode collect_variable_data(std::vector<std::string>& var_names, std::vector<int>& tstep_nums,
+ std::vector<double>& tstep_vals, EntityHandle fileSet);
- // initialize file: this is where all defines are done
+ // Initialize file: this is where all defines are done
// the VarData dimension ids are filled up after define
- ErrorCode initialize_file( std::vector<std::string> & var_names); // these are from options
+ ErrorCode initialize_file(std::vector<std::string>& var_names); // These are from options
- // take the info from VarData and write first the coordinates, then the actual variables
- ErrorCode write_values(std::vector<std::string> & var_names, EntityHandle fileSet);
- // interface instance
- Interface *mbImpl;
+ // Take the info from VarData and write first the coordinates, then the actual variables
+ ErrorCode write_values(std::vector<std::string>& var_names, EntityHandle fileSet);
+
+ template <typename T> void jik_to_kji(size_t ni, size_t nj, size_t nk, T* dest, T* source)
+ {
+ size_t nik = ni * nk, nij = ni * nj;
+ for (std::size_t k = 0; k != nk; k++)
+ for (std::size_t j = 0; j != nj; j++)
+ for (std::size_t i = 0; i != ni; i++)
+ dest[k*nij + j*ni + i] = source[j*nik + i*nk + k];
+ }
+
+ // Interface instance
+ Interface* mbImpl;
WriteUtilIface* mWriteIface;
- // File var
- const char *fileName;
+ //! File var
+ const char* fileName;
int IndexFile;
//! File numbers assigned by (p)netcdf
int fileId;
-
//! Debug stuff
DebugOutput dbgOut;
@@ -200,29 +199,29 @@ private:
#ifdef USE_MPI
ParallelComm* myPcomm;
#endif
- //! write options
- //
+
+ //! Write options
bool noMesh;
bool noVars;
+
/*
- * not used yet, maybe later
+ * Not used yet, maybe later
bool spectralMesh;
bool noMixedElements;
bool noEdges;*/
int gatherSetRank;
//! Cached tags for writing. this will be important for ordering the data, in parallel
- //
Tag mGlobalIdTag;
//! Are we writing in parallel? (probably in the future)
bool isParallel;
- // CAM Euler, etc,
+ //! CAM Euler, etc,
std::string grid_type;
+
//! Helper class instance
NCWriteHelper* myHelper;
-
};
} // namespace moab
Repository URL: https://bitbucket.org/fathomteam/moab/
--
This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.
More information about the moab-dev
mailing list