From 924f0a6d40bd4fc09715e37f8f2029fcef95c908 Mon Sep 17 00:00:00 2001 From: Thorsten Liebig Date: Mon, 5 Apr 2010 20:22:03 +0200 Subject: [PATCH] first hdf5-file dump implementation --- FDTD/processfields.cpp | 116 ++++++++++++++++++++++++++++++++++-- FDTD/processfields.h | 28 +++++++-- FDTD/processfields_td.cpp | 105 +++++++++++++++++++++++++------- FDTD/processfields_td.h | 4 +- matlab/examples/PlaneWave.m | 8 +-- openEMS.pro | 4 +- openems.cpp | 7 ++- 7 files changed, 232 insertions(+), 40 deletions(-) diff --git a/FDTD/processfields.cpp b/FDTD/processfields.cpp index ff70810..e988ff6 100644 --- a/FDTD/processfields.cpp +++ b/FDTD/processfields.cpp @@ -17,10 +17,14 @@ #include "processfields.h" +#include "H5Cpp.h" + ProcessFields::ProcessFields(Operator* op, Engine* eng) : Processing(op, eng) { - DumpMode=0; - DumpType = 0; + m_DumpMode = NO_INTERPOLATION; + m_DumpType = E_FIELD_DUMP; + // vtk-file is default + m_fileType = VTK_FILETYPE; // SetSubSampling(1); for (int n=0;n<3;++n) @@ -42,9 +46,68 @@ ProcessFields::~ProcessFields() } } +void ProcessFields::InitProcess() +{ + string names[] = {"x","y","z"}; + if (m_fileType==HDF5_FILETYPE) + { + unsigned int* NrLines; + double** Lines; + + if (m_DumpMode==CELL_INTERPOLATE) + { + NrLines = numDLines; + Lines = discDLines; + } + else if (m_DumpMode==NO_INTERPOLATION) + { + NrLines = numLines; + Lines = discLines; + } + else + return; + + m_fileName+= ".h5"; + H5::H5File* file = new H5::H5File( m_fileName , H5F_ACC_TRUNC ); + + H5::Group* group = new H5::Group( file->createGroup( "/Mesh" )); + for (int n=0;n<3;++n) + { + hsize_t dimsf[1]; // dataset dimensions + dimsf[0] = NrLines[n]; + H5::DataSpace dataspace( 1, dimsf ); + H5::FloatType datatype( H5::PredType::NATIVE_FLOAT ); + H5::DataSet dataset = group->createDataSet( names[n].c_str(), datatype, dataspace ); + //convert to float... + float* array = new float[NrLines[n]]; + for (int i=0;icreateGroup( "/FieldData" )); + delete group; + delete file; + } +} + +string ProcessFields::GetFieldNameByType(DumpType type) +{ + switch (type) + { + case E_FIELD_DUMP: + return "E-Field"; + case H_FIELD_DUMP: + return "H-Field"; + } + return "unknown field"; +} + void ProcessFields::DefineStartStopCoord(double* dstart, double* dstop) { - if (DumpMode==0) + if (m_DumpMode==NO_INTERPOLATION) { if (Op->SnapToMesh(dstart,start)==false) cerr << "ProcessFields::DefineStartStopCoord: Warning: Snapping problem, check start value!!" << endl; if (Op->SnapToMesh(dstop,stop)==false) cerr << "ProcessFields::DefineStartStopCoord: Warning: Snapping problem, check stop value!!" << endl; @@ -69,7 +132,7 @@ void ProcessFields::DefineStartStopCoord(double* dstart, double* dstop) } } } - else if (DumpMode==2) + else if (m_DumpMode==CELL_INTERPOLATE) { if (Op->SnapToMesh(dstart,start,true)==false) cerr << "ProcessFields::DefineStartStopCoord: Warning: Snapping problem, check start value!!" << endl; if (Op->SnapToMesh(dstop,stop,true)==false) cerr << "ProcessFields::DefineStartStopCoord: Warning: Snapping problem, check stop value!!" << endl; @@ -86,7 +149,7 @@ void ProcessFields::DefineStartStopCoord(double* dstart, double* dstop) } ++stop[n]; numDLines[n]=stop[n]-start[n]; - // cerr << " number of lines " << numDLines[n] << endl; +// cerr << " number of lines " << numDLines[n] << endl; delete[] discDLines[n]; discDLines[n] = new double[numDLines[n]]; for (unsigned int i=0;i "tmp/efield_000045.vtk" for timestep 45 or "tmp/efield_2.40000e9.vtk" for 2.4GHz E-field dump. + //! Used file pattern e.g. pattern="tmp/efield_" --> "tmp/efield_000045.vtk" for timestep 45 or "tmp/efield_2.40000e9.vtk" for 2.4GHz E-field dump. (VTK FileType only) \sa SetFileType void SetFilePattern(string fp) {filePattern=fp;} + //! Set the filename for a hdf5 data group file (HDF5 FileType only) \sa SetFileType + void SetFileName(string fn) {m_fileName=fn;} + //! Define the Dump-Mode - void SetDumpMode(int mode) {DumpMode=mode;} + void SetDumpMode(DumpMode mode) {m_DumpMode=mode;} //! This methode will dump all fields in the center of a main cell (dual-node) using 4 E-field and 2 H-fields per direction. (default) - void SetDumpMode2Cell() {DumpMode=2;} + void SetDumpMode2Cell() {m_DumpMode=CELL_INTERPOLATE;} //! Set dump type: 0 for E-fields, 1 for H-fields, 2 for D-fields, 3 for B-fields, 4 for J-fields, etc... - void SetDumpType(int type) {DumpType=type;} + void SetDumpType(DumpType type) {m_DumpType=type;} static bool DumpVectorArray2VTK(ofstream &file, string name, FDTD_FLOAT**** array, double** discLines, unsigned int* numLines); static bool DumpMultiVectorArray2VTK(ofstream &file, string names[], FDTD_FLOAT**** array[], unsigned int numFields, double** discLines, unsigned int* numLines); static bool DumpScalarArray2VTK(ofstream &file, string name, FDTD_FLOAT*** array, double** discLines, unsigned int* numLines); static bool DumpMultiScalarArray2VTK(ofstream &file, string names[], FDTD_FLOAT*** array[], unsigned int numFields, double** discLines, unsigned int* numLines); + static bool DumpVectorArray2HDF5(string filename, string name, FDTD_FLOAT**** array, unsigned int* numLines); + double CalcTotalEnergy(); + void SetFileType(FileType fileType) {m_fileType=fileType;} + // virtual void Process(); protected: static void WriteVTKHeader(ofstream &file, double** discLines, unsigned int* numLines); static void WriteVTKVectorArray(ofstream &file, string name, FDTD_FLOAT**** array, unsigned int* numLines); static void WriteVTKScalarArray(ofstream &file, string name, FDTD_FLOAT*** array, unsigned int* numLines); - int DumpMode; - int DumpType; + static string GetFieldNameByType(DumpType type); + DumpMode m_DumpMode; + DumpType m_DumpType; string filePattern; + string m_fileName; + FileType m_fileType; // unsigned int subSample[3]; diff --git a/FDTD/processfields_td.cpp b/FDTD/processfields_td.cpp index bc8a93c..25b476a 100644 --- a/FDTD/processfields_td.cpp +++ b/FDTD/processfields_td.cpp @@ -29,12 +29,12 @@ ProcessFieldsTD::~ProcessFieldsTD() { } -void ProcessFieldsTD::DumpCellInterpol(ofstream &file) +void ProcessFieldsTD::DumpCellInterpol(string filename) { FDTD_FLOAT**** volt = Eng->GetVoltages(); FDTD_FLOAT**** curr = Eng->GetCurrents(); - if (DumpType==0) + if (m_DumpType==E_FIELD_DUMP) { //create array FDTD_FLOAT**** E_T = Create_N_3DArray(numDLines); @@ -66,12 +66,27 @@ void ProcessFieldsTD::DumpCellInterpol(ofstream &file) } } } - DumpVectorArray2VTK(file,string("E-Field"),E_T,discDLines,numDLines); + + if (m_fileType==VTK_FILETYPE) + { + ofstream file(filename.c_str()); + if (file.is_open()==false) { cerr << "ProcessFieldsTD::Process: can't open file '" << filename << "' for writing... abort! " << endl;}; + DumpVectorArray2VTK(file,string("E-Field"),E_T,discDLines,numDLines); + file.close(); + } + else if (m_fileType==HDF5_FILETYPE) + { + stringstream ss; + ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); + DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),E_T,numDLines); + } + else + cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl; Delete_N_3DArray(E_T,numDLines); E_T = NULL; } - if (DumpType==1) + if (m_DumpType==1) { //create array FDTD_FLOAT**** H_T = Create_N_3DArray(numDLines); @@ -104,20 +119,34 @@ void ProcessFieldsTD::DumpCellInterpol(ofstream &file) } } } - DumpVectorArray2VTK(file,string("H-Field"),H_T,discDLines,numDLines); + if (m_fileType==VTK_FILETYPE) + { + ofstream file(filename.c_str()); + if (file.is_open()==false) { cerr << "ProcessFieldsTD::Process: can't open file '" << filename << "' for writing... abort! " << endl;}; + DumpVectorArray2VTK(file,string("H-Field"),H_T,discDLines,numDLines); + file.close(); + } + else if (m_fileType==HDF5_FILETYPE) + { + stringstream ss; + ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); + DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),H_T,numDLines); + } + else + cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl; Delete_N_3DArray(H_T,numDLines); H_T = NULL; } } -void ProcessFieldsTD::DumpNoInterpol(ofstream &file) +void ProcessFieldsTD::DumpNoInterpol(string filename) { FDTD_FLOAT**** volt = Eng->GetVoltages(); FDTD_FLOAT**** curr = Eng->GetCurrents(); unsigned int pos[3]; double delta[3]; - if (DumpType==0) + if (m_DumpType==E_FIELD_DUMP) { //create array FDTD_FLOAT**** E_T = Create_N_3DArray(numLines); @@ -136,12 +165,27 @@ void ProcessFieldsTD::DumpNoInterpol(ofstream &file) } } } - DumpVectorArray2VTK(file,string("E-Field"),E_T,discLines,numLines); + if (m_fileType==VTK_FILETYPE) + { + ofstream file(filename.c_str()); + if (file.is_open()==false) { cerr << "ProcessFieldsTD::Process: can't open file '" << filename << "' for writing... abort! " << endl;}; + DumpVectorArray2VTK(file,string("E-Field"),E_T,discLines,numLines); + file.close(); + } + else if (m_fileType==HDF5_FILETYPE) + { + stringstream ss; + ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); + DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),E_T,numLines); + } + else + cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl; + Delete_N_3DArray(E_T,numLines); E_T = NULL; } - if (DumpType==1) + if (m_DumpType==H_FIELD_DUMP) { //create array FDTD_FLOAT**** H_T = Create_N_3DArray(numLines); @@ -161,7 +205,22 @@ void ProcessFieldsTD::DumpNoInterpol(ofstream &file) } } } - DumpVectorArray2VTK(file,string("H-Field"),H_T,discLines,numLines); + if (m_fileType==VTK_FILETYPE) + { + ofstream file(filename.c_str()); + if (file.is_open()==false) { cerr << "ProcessFieldsTD::Process: can't open file '" << filename << "' for writing... abort! " << endl;}; + DumpVectorArray2VTK(file,string("H-Field"),H_T,discLines,numLines); + file.close(); + } + else if (m_fileType==HDF5_FILETYPE) + { + stringstream ss; + ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); + DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),H_T,numLines); + } + else + cerr << "ProcessFieldsTD::DumpCellInterpol: unknown File-Type" << endl; + Delete_N_3DArray(H_T,numLines); H_T = NULL; } @@ -173,16 +232,22 @@ int ProcessFieldsTD::Process() if (filePattern.empty()) return -1; if (CheckTimestep()==false) return GetNextInterval(); stringstream ss; - ss << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); + ss << filePattern << std::setw( pad_length ) << std::setfill( '0' ) << Eng->GetNumberOfTimesteps(); - string filename = filePattern + ss.str() + ".vtk"; - ofstream file(filename.c_str()); - if (file.is_open()==false) { cerr << "ProcessFieldsTD::Process: can't open file '" << filename << "' for writing... abort! " << endl; return GetNextInterval();}; - - if (DumpMode==0) - DumpNoInterpol(file); - if (DumpMode==2) - DumpCellInterpol(file); - file.close(); + if (m_fileType==VTK_FILETYPE) + { + ss << ".vtk"; + if (m_DumpMode==NO_INTERPOLATION) + DumpNoInterpol(ss.str()); + if (m_DumpMode==CELL_INTERPOLATE) + DumpCellInterpol(ss.str()); + } + else if (m_fileType==HDF5_FILETYPE) + { + if (m_DumpMode==NO_INTERPOLATION) + DumpNoInterpol(m_fileName); + if (m_DumpMode==CELL_INTERPOLATE) + DumpCellInterpol(m_fileName); + } return GetNextInterval(); } diff --git a/FDTD/processfields_td.h b/FDTD/processfields_td.h index 52ea869..84bfc63 100644 --- a/FDTD/processfields_td.h +++ b/FDTD/processfields_td.h @@ -34,8 +34,8 @@ public: protected: int pad_length; - void DumpNoInterpol(ofstream &file); - void DumpCellInterpol(ofstream &file); + void DumpNoInterpol(string filename); + void DumpCellInterpol(string filename); }; #endif // PROCESSFIELDS_TD_H diff --git a/matlab/examples/PlaneWave.m b/matlab/examples/PlaneWave.m index 25e531e..3717ec7 100644 --- a/matlab/examples/PlaneWave.m +++ b/matlab/examples/PlaneWave.m @@ -53,13 +53,13 @@ CSX = AddExcitation(CSX,'excite',0,[0 1 0]); CSX = AddBox(CSX,'excite',0 ,start,stop); %dump -CSX = AddDump(CSX,'Et_',0,0); +CSX = AddDump(CSX,'Et',0,0,1); start = [mesh.x(1) , 0 , mesh.z(1)]; stop = [mesh.x(end) , 0 , mesh.z(end)]; -CSX = AddBox(CSX,'Et_',0 , start,stop); +CSX = AddBox(CSX,'Et',0 , start,stop); -CSX = AddDump(CSX,'Ht_',1,0); -CSX = AddBox(CSX,'Ht_',0,start,stop); +CSX = AddDump(CSX,'Ht',1,0,1); +CSX = AddBox(CSX,'Ht',0,start,stop); %Write openEMS compatoble xml-file WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX); diff --git a/openEMS.pro b/openEMS.pro index 34c098e..d66d34d 100644 --- a/openEMS.pro +++ b/openEMS.pro @@ -15,10 +15,12 @@ LIBS += -L../CSXCAD \ -lfparser \ -L../tinyxml \ -ltinyxml \ - -lboost_thread + -lboost_thread \ + -lhdf5_cpp QMAKE_LFLAGS += \'-Wl,-rpath,\$$ORIGIN/../CSXCAD\' QMAKE_LFLAGS += \'-Wl,-rpath,\$$ORIGIN/../fparser\' QMAKE_LFLAGS += \'-Wl,-rpath,\$$ORIGIN/../tinyxml\' + SOURCES += main.cpp \ tools/ErrorMsg.cpp \ tools/AdrOp.cpp \ diff --git a/openems.cpp b/openems.cpp index 79f06ae..f41360c 100644 --- a/openems.cpp +++ b/openems.cpp @@ -346,10 +346,13 @@ int openEMS::SetupFDTD(const char* file) CSPropDumpBox* db = DumpProps.at(i)->ToDumpBox(); if (db) { - ProcTD->SetDumpType(db->GetDumpType()); - ProcTD->SetDumpMode(db->GetDumpMode()); + ProcTD->SetDumpType((ProcessFields::DumpType)db->GetDumpType()); + ProcTD->SetDumpMode((ProcessFields::DumpMode)db->GetDumpMode()); + ProcTD->SetFileType((ProcessFields::FileType)db->GetFileType()); ProcTD->SetFilePattern(db->GetName()); + ProcTD->SetFileName(db->GetName()); ProcTD->DefineStartStopCoord(start,stop); + ProcTD->InitProcess(); PA->AddProcessing(ProcTD); } else