Skip to content

Commit

Permalink
add in VTK header fields 6 1D arrays with native coordinates at edges…
Browse files Browse the repository at this point in the history
… and centers
  • Loading branch information
volodia99 committed Nov 5, 2024
1 parent d6960c9 commit 9f25ea3
Show file tree
Hide file tree
Showing 3 changed files with 91 additions and 70 deletions.
49 changes: 34 additions & 15 deletions pytools/vtk_io.py
Original file line number Diff line number Diff line change
Expand Up @@ -60,6 +60,7 @@ def _load_header(self, fh, geometry=None):

self.geometry = geometry
self._has_native_coordinates = False
self.native_coords = {}

ref_position = fh.tell()
line = fh.readline().decode("utf-8") # DIMENSIONS NX NY NZ or FIELD
Expand Down Expand Up @@ -89,8 +90,25 @@ def _load_header(self, fh, geometry=None):
self.t = np.fromfile(fh, dt, 1)
elif d.startswith("PERIODICITY"):
self.periodicity = np.fromfile(fh, dtype=dint, count=3).astype(bool)
elif d.startswith("HAS_NATIVE_COORDINATES"):
self._has_native_coordinates = np.fromfile(fh, dtype=dint, count=1).astype(bool)
elif d.startswith("X1L_NATIVE_COORDINATES"):
self._has_native_coordinates = True
native_name, native_dim = (d.split()[0], int(d.split()[2]))
self.native_coords[native_name] = np.fromfile(fh, dtype=dt, count=native_dim)
elif d.startswith("X2L_NATIVE_COORDINATES"):
native_name, native_dim = (d.split()[0], int(d.split()[2]))
self.native_coords[native_name] = np.fromfile(fh, dtype=dt, count=native_dim)
elif d.startswith("X3L_NATIVE_COORDINATES"):
native_name, native_dim = (d.split()[0], int(d.split()[2]))
self.native_coords[native_name] = np.fromfile(fh, dtype=dt, count=native_dim)
elif d.startswith("X1C_NATIVE_COORDINATES"):
native_name, native_dim = (d.split()[0], int(d.split()[2]))
self.native_coords[native_name] = np.fromfile(fh, dtype=dt, count=native_dim)
elif d.startswith("X2C_NATIVE_COORDINATES"):
native_name, native_dim = (d.split()[0], int(d.split()[2]))
self.native_coords[native_name] = np.fromfile(fh, dtype=dt, count=native_dim)
elif d.startswith("X3C_NATIVE_COORDINATES"):
native_name, native_dim = (d.split()[0], int(d.split()[2]))
self.native_coords[native_name] = np.fromfile(fh, dtype=dt, count=native_dim)
else:
warnings.warn("Found unknown field %s" % d)
fh.readline() # skip extra linefeed (empty line)
Expand Down Expand Up @@ -344,25 +362,26 @@ def _load_hydro(self, fh):
if self._has_native_coordinates:
if self.geometry == "spherical":
native2attr = {
"X1L_NATIVE_COORDS": "rl",
"X1C_NATIVE_COORDS": "r",
"X2L_NATIVE_COORDS": "thetal",
"X2C_NATIVE_COORDS": "theta",
"X3L_NATIVE_COORDS": "phil",
"X3C_NATIVE_COORDS": "phi",
"X1L_NATIVE_COORDINATES": "rl",
"X1C_NATIVE_COORDINATES": "r",
"X2L_NATIVE_COORDINATES": "thetal",
"X2C_NATIVE_COORDINATES": "theta",
"X3L_NATIVE_COORDINATES": "phil",
"X3C_NATIVE_COORDINATES": "phi",
}
elif self.geometry in ("cartesian", "cylindrical", "polar"):
native2attr = {
"X1L_NATIVE_COORDS": "xl",
"X1C_NATIVE_COORDS": "x",
"X2L_NATIVE_COORDS": "yl",
"X2C_NATIVE_COORDS": "y",
"X3L_NATIVE_COORDS": "zl",
"X3C_NATIVE_COORDS": "z",
"X1L_NATIVE_COORDINATES": "xl",
"X1C_NATIVE_COORDINATES": "x",
"X2L_NATIVE_COORDINATES": "yl",
"X2C_NATIVE_COORDINATES": "y",
"X3L_NATIVE_COORDINATES": "zl",
"X3C_NATIVE_COORDINATES": "z",
}

for native_field, attr in native2attr.items():
setattr(self, attr, self.data[native_field])
np.testing.assert_array_equal(getattr(self, attr), self.native_coords[native_field])
setattr(self, attr, self.native_coords[native_field])

def _load_particles(self, fh):
raise NotImplementedError("Particles vtk are not supported yet !")
Expand Down
111 changes: 57 additions & 54 deletions src/output/vtk.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -330,8 +330,6 @@ int Vtk::Write() {
WriteScalar(fileHdl, vect3D, name);
}

WriteNativeCoords(fileHdl);

#ifdef WITH_MPI
MPI_SAFE_CALL(MPI_File_close(&fileHdl));
#else
Expand Down Expand Up @@ -389,8 +387,8 @@ void Vtk::WriteHeader(IdfxFileHandler fvtk, real time) {
#elif VTK_FORMAT == VTK_STRUCTURED_GRID
ssheader << "DATASET STRUCTURED_GRID" << std::endl;
#endif
// fields: geometry, periodicity, time, hasNativeCoordinates
int nfields = 4;
// fields: geometry, periodicity, time, 6 NativeCoordinates (x1l, x2l, x3l, x1, x2, x3)
int nfields = 9;

// Write grid geometry in the VTK file
ssheader << "FIELD FieldData " << nfields << std::endl;
Expand Down Expand Up @@ -440,18 +438,69 @@ void Vtk::WriteHeader(IdfxFileHandler fvtk, real time) {
// Done, add cariage return for next ascii write
ssheader << std::endl;

// write bool if native coordinates
ssheader << "HAS_NATIVE_COORDINATES 1 1 int" << std::endl;
// write x1l native coordinates
ssheader << "X1L_NATIVE_COORDINATES 1 " << (nx1 + ioffset) << " float" << std::endl;
// Flush the ascii header
header = ssheader.str();
WriteHeaderString(header.c_str(), fvtk);
// reset the string stream
ssheader.str(std::string());
WriteHeaderBinary(this->xnode, nx1 + ioffset, fvtk);
// Done, add cariage return for next ascii write
ssheader << std::endl;

// write x2l native coordinates
ssheader << "X2L_NATIVE_COORDINATES 1 " << (nx2 + joffset) << " float" << std::endl;
// Flush the ascii header
header = ssheader.str();
WriteHeaderString(header.c_str(), fvtk);
// reset the string stream
ssheader.str(std::string());
WriteHeaderBinary(this->ynode, nx2 + joffset, fvtk);
// Done, add cariage return for next ascii write
ssheader << std::endl;

// write x3l native coordinates
ssheader << "X3L_NATIVE_COORDINATES 1 " << (nx3 + koffset) << " float" << std::endl;
// Flush the ascii header
header = ssheader.str();
WriteHeaderString(header.c_str(), fvtk);
// reset the string stream
ssheader.str(std::string());
WriteHeaderBinary(this->znode, nx3 + koffset, fvtk);
// Done, add cariage return for next ascii write
ssheader << std::endl;

int32_t hasNativeCoordinatesBig = bigEndian(1);
WriteHeaderBinary(&hasNativeCoordinatesBig, 1, fvtk);
// write x1 native coordinates
ssheader << "X1C_NATIVE_COORDINATES 1 " << nx1 << " float" << std::endl;
// Flush the ascii header
header = ssheader.str();
WriteHeaderString(header.c_str(), fvtk);
// reset the string stream
ssheader.str(std::string());
WriteHeaderBinary(this->xcenter, nx1, fvtk);
// Done, add cariage return for next ascii write
ssheader << std::endl;

// write x2 native coordinates
ssheader << "X2C_NATIVE_COORDINATES 1 " << nx2 << " float" << std::endl;
// Flush the ascii header
header = ssheader.str();
WriteHeaderString(header.c_str(), fvtk);
// reset the string stream
ssheader.str(std::string());
WriteHeaderBinary(this->ycenter, nx2, fvtk);
// Done, add cariage return for next ascii write
ssheader << std::endl;

// write x3 native coordinates
ssheader << "X3C_NATIVE_COORDINATES 1 " << nx3 << " float" << std::endl;
// Flush the ascii header
header = ssheader.str();
WriteHeaderString(header.c_str(), fvtk);
// reset the string stream
ssheader.str(std::string());
WriteHeaderBinary(this->zcenter, nx3, fvtk);
// Done, add cariage return for next ascii write
ssheader << std::endl;

Expand Down Expand Up @@ -513,52 +562,6 @@ void Vtk::WriteHeader(IdfxFileHandler fvtk, real time) {
#undef VTK_STRUCTURED_GRID
#undef VTK_RECTILINEAR_GRID

/* ********************************************************************* */
void Vtk::WriteNativeCoords(IdfxFileHandler fvtk) {
/*!
* Write native coordinates in VTK.
* \param [in] fvtk pointer to file
*/
std::stringstream native_coordx1l, native_coordx1;
std::stringstream native_coordx2l, native_coordx2;
std::stringstream native_coordx3l, native_coordx3;
std::string header;

native_coordx1l << "X1L_NATIVE_COORDS " << this->nx1 + this->ioffset;
native_coordx1l << " float" << std::endl;
header = native_coordx1l.str();
WriteHeaderString(header.c_str(), fvtk);
WriteHeaderBinary(xnode, nx1 + ioffset, fvtk);

native_coordx1 << "X1C_NATIVE_COORDS " << this->nx1;
native_coordx1 << " float" << std::endl;
header = native_coordx1.str();
WriteHeaderString(header.c_str(), fvtk);
WriteHeaderBinary(xcenter, nx1, fvtk);

native_coordx2l << "X2L_NATIVE_COORDS " << this->nx2 + this->joffset;
native_coordx2l << " float" << std::endl;
header = native_coordx2l.str();
WriteHeaderString(header.c_str(), fvtk);
WriteHeaderBinary(ynode, nx2 + joffset, fvtk);

native_coordx2 << "X2C_NATIVE_COORDS " << this->nx2 << " float" << std::endl;
header = native_coordx2.str();
WriteHeaderString(header.c_str(), fvtk);
WriteHeaderBinary(ycenter, nx2, fvtk);

native_coordx3l << "X3L_NATIVE_COORDS " << this->nx3 + this->koffset;
native_coordx3l << " float" << std::endl;
header = native_coordx3l.str();
WriteHeaderString(header.c_str(), fvtk);
WriteHeaderBinary(znode, nx3 + koffset, fvtk);

native_coordx3 << "X3C_NATIVE_COORDS " << this->nx3 << " float" << std::endl;
header = native_coordx3.str();
WriteHeaderString(header.c_str(), fvtk);
WriteHeaderBinary(zcenter, nx3, fvtk);
}

/* ********************************************************************* */
void Vtk::WriteScalar(IdfxFileHandler fvtk, float* Vin, const std::string &var_name) {
/*!
Expand Down
1 change: 0 additions & 1 deletion src/output/vtk.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -149,7 +149,6 @@ class Vtk : public BaseVtk {
void WriteHeader(IdfxFileHandler, real);
void WriteScalar(IdfxFileHandler, float*, const std::string &);
void WriteHeaderNodes(IdfxFileHandler);
void WriteNativeCoords(IdfxFileHandler);

// output directory
fs::path outputDirectory;
Expand Down

0 comments on commit 9f25ea3

Please sign in to comment.