User: reisswig Date: 2013/06/04 03:40 PM
Modified: /PITTNullCode/SphericalHarmonicReconGen/ param.ccl /PITTNullCode/SphericalHarmonicReconGen/src/ sph_database.cc, sph_database.hh, startup.cc
Log: SphericalHarmonicReconGen: Now supports SpEC hdf5-v2 format
File Changes:
Directory: /PITTNullCode/SphericalHarmonicReconGen/ ===================================================
File [modified]: param.ccl Delta lines: +4 -3 =================================================================== --- PITTNullCode/SphericalHarmonicReconGen/param.ccl 2013-05-27 17:02:26 UTC (rev 84) +++ PITTNullCode/SphericalHarmonicReconGen/param.ccl 2013-06-04 20:40:27 UTC (rev 85) @@ -5,9 +5,10 @@
KEYWORD format "Boundary data file format" { - ASCII :: "standard ASCII format: lm-modes and spheres are stored in a 2d array (gnuplot style) for a given timestep" - DAT :: "optimized ASCII format: all variables and lm-modes are stored in one row for a given timestep" - SpEC-H5 :: "SpEC-HDF5 format: Data is stored in HDF5 files grouped according to SpEC output." + ASCII :: "standard ASCII format: lm-modes and spheres are stored in a 2d array (gnuplot style) for a given timestep" + DAT :: "optimized ASCII format: all variables and lm-modes are stored in one row for a given timestep" + SpEC-H5 :: "SpEC-HDF5 format: Data is stored in HDF5 files grouped according to SpEC output." + SpEC-H5-v2 :: "SpEC-HDF5 format (v2): Data is stored in HDF5 files grouped according to new SpEC output." } "ASCII"
CCTK_INT sphere_number "the sphere number (if multiple spheres are present in one file) that corresponds to the worldtube radius NullSHRExtract::cr"
Directory: /PITTNullCode/SphericalHarmonicReconGen/src/ =======================================================
File [modified]: sph_database.cc Delta lines: +224 -0 =================================================================== --- PITTNullCode/SphericalHarmonicReconGen/src/sph_database.cc 2013-05-27 17:02:26 UTC (rev 84) +++ PITTNullCode/SphericalHarmonicReconGen/src/sph_database.cc 2013-06-04 20:40:27 UTC (rev 85) @@ -1066,7 +1066,231 @@
+//------------------------------------------------------------------------------------------- +// +// +// SpEC-H5-v2 format goes here........... +// +// +//------------------------------------------------------------------------------------------
+ +/*const string SPH_db_SpEC_H5_v2::time_attrib_name = "Time"; +const string SPH_db_SpEC_H5_v2::step_basename = "Step"; +*/ + // table of varibale names: this maps a column number to a group. + // The first entry per varibale is the main group, the second is a possible subgroup, e.g. + // /g/Step00001/xx + // for SpEC +const string SPH_db_SpEC_H5_v2::varname_table[30] = + { "Lapse.dat", + "Shiftx.dat" , "Shifty.dat", "Shiftz.dat", + "gxx.dat", "gxy.dat", "gxz.dat", + "gyy.dat", "gyz.dat", + "gzz.dat", + + "DrLapse.dat", + "DrShiftx.dat", "DrShifty.dat", "DrShiftz.dat", + "Drgxx.dat", "Drgxy.dat", "Drgxz.dat", + "Drgyy.dat", "Drgyz.dat", + "Drgzz.dat", + + "DtLapse.dat", + "DtShiftx.dat", "DtShifty.dat" "DtShiftz.dat", + "Dtgxx.dat", "Dtgxy.dat", "Dtgxz.dat", + "Dtgyy.dat", "Dtgyz.dat", + "Dtgzz.dat" }; + + + +SPH_db_SpEC_H5_v2::SPH_db_SpEC_H5_v2(const string& fname_, + const int cached_timesteps) + : SPH_database(cached_timesteps), _fname(fname_) +{ + // open file + file = H5Fopen(_fname.c_str(), H5F_ACC_RDONLY, H5P_DEFAULT); + + scan(); } + + +SPH_db_SpEC_H5_v2::~SPH_db_SpEC_H5_v2() +{ + // close file + H5Fclose(file); +} + + +int SPH_db_SpEC_H5_v2::remove_non_monotonic_steps(const int k) +{ + // go backwards from k and remove until we get monotonicity + const double time = steps[k]; + const int size = steps.size(); + int c = 1; + for (int i=k-1; i >= 0; --i) { + if (time > steps[i]) { + break; + } + ++c; + } + + steps.erase(steps.begin()+k-c, steps.begin()+k); + _times.erase(_times.begin()+k-c, _times.begin()+k); + _iterations.erase(_iterations.begin()+k-c, _iterations.begin()+k); + + // return the number of steps that were removed + return size-steps.size(); +} + + + +void SPH_db_SpEC_H5_v2::scan() +{ + + // scan datasets + scan_HDF5(); + + + // check continuity and delta_t of timesteps + if (_times.size() > 1) + { + cout << "SPH_db_SpEC_H5 " << _fname << ": t_0 = " << _times[0] << ", t_final = " << _times.back() << endl; + _delta_t = 0; + for (size_t i=1; i < _times.size(); ++i) + { + if (fabs(_times[i] - _times[i-1]) < 1e-8 || _times[i] < _times[i-1]) { + ostringstream str; + str << "Removing non-monotonic timesteps: time[" << i-1 << "] = " << _times[i-1] << ", time[" << i << "] = " << _times[i]; + CCTK_WARN(1, str.str().c_str()); + i -= remove_non_monotonic_steps(i); + continue; + } + if (fabs(_times[i] - _times[i-1] - _delta_t) >= 1e-8) + { + _delta_t = _times[i]-_times[i-1]; + cout << "SPH_db_SpEC_H5 " << _fname << ": t = " << _times[i-1] << ", delta_t = " << _delta_t << endl; + } + } + // set _delta_t to correspond to initial delta_t + _delta_t = _times[1]-_times[0]; + } + +} + + +void SPH_db_SpEC_H5_v2::read(const size_t timestep, const size_t varno, + vector<vector<vector<CCTK_COMPLEX> > >& coeff) const +{ + if (timestep >= steps.size()) { + CCTK_WARN(0, "Requested timestep not in worldtube data file. Stopping!"); + } + + hid_t dataset = H5Dopen(file, varname_table[varno].c_str());//, H5P_DEFAULT); + vector<double> buffer((_lmax+1)*(_lmax+1)*2); + + hsize_t dimsm[1] = { buffer.size() }; + hid_t memspace = H5Screate_simple (1, dimsm, NULL); + + hid_t dataspace = H5Dget_space(dataset); + hsize_t dims[2]; + H5Sget_simple_extent_dims (dataspace, dims, NULL); + + hsize_t offset[2] = { timestep, 1 }; + hsize_t count[2] = { 1, buffer.size() }; + H5Sselect_hyperslab (dataspace, H5S_SELECT_SET, offset, NULL, + count, NULL); + + H5Dread(dataset, H5T_NATIVE_DOUBLE, memspace, dataspace, H5P_DEFAULT, &buffer.front()); + + H5Sclose(memspace); + H5Sclose(dataspace); + H5Dclose(dataset); + + int c = 0; + for (int l=0; l <= _lmax; ++l) { + for (int m=l; m >= -l; --m, ++c) { + //coeff[0][l][m+l] = CCTK_Cmplx(buffer[c], buffer[c + buffer.size()/2]); + coeff[0][l][m+l] = CCTK_Cmplx(buffer[2*c], buffer[2*c + 1]); + } + } +} + + + +void SPH_db_SpEC_H5_v2::scan_HDF5() +{/* + // pick one variable and browse through entire file + // to find out how many (and what) cycles and times there are + // (assuming that all variables are present at the same timesteps) + H5Literate_by_name (file, varname_table[0][0].c_str(), H5_INDEX_NAME, H5_ITER_NATIVE, NULL, H5iter, this, H5P_DEFAULT); + + // sort 'steps' vector + //sort(steps.begin(), steps.end()); + _iterations.resize(steps.size()); + _times.resize(steps.size()); + for (size_t i=0; i < steps.size(); ++i) { + //cout << steps[i].name << " " << steps[i].time << endl; + _iterations[i] = i; + _times[i] = steps[i].time; + } +*/ + + // get lmax from first variable and timestep + hid_t dataset = H5Dopen(file, "gxx.dat"); //, H5P_DEFAULT); + //cout << (varname_table[0][0]+"/"+steps[0].name+"/"+varname_table[0][1]) << endl; + hid_t dataspace = H5Dget_space(dataset); + hsize_t dims[2]; + H5Sget_simple_extent_dims (dataspace, dims, NULL); + + // first dim is number of timesteps + steps.resize(dims[0]); + // second dim is the number of lmodes + _lmax = sqrt(dims[1]/2)-1; + + // select hyperslab in dataset + hsize_t offset[2] = { 0, 0 }; + hsize_t count[2] = { dims[0], 1 }; + H5Sselect_hyperslab (dataspace, H5S_SELECT_SET, offset, NULL, + count, NULL); + + // define dataspace of memory + hsize_t dimsm[1] = { steps.size() }; + hid_t memspace = H5Screate_simple (1, dimsm, NULL); + + // read "times" + H5Dread (dataset, H5T_NATIVE_DOUBLE, memspace, dataspace, + H5P_DEFAULT, &steps.front()); + + H5Sclose(memspace); + H5Sclose(dataspace); + H5Dclose(dataset); + + _n_variables = 30; + _n_spheres = 1; + + _iterations.resize(steps.size()); + _times.resize(steps.size()); + for (size_t i=0; i < steps.size(); ++i) { + //cout << steps[i].name << " " << steps[i].time << endl; + _iterations[i] = i; + _times[i] = steps[i]; + } + + // cross-check that all expected variables are present + + + // display some info + cout << "SPH_db_SpEC_H5_v2 " << _fname << ": n_spheres = " << _n_spheres << endl; + cout << "SPH_db_SpEC_H5_v2 " << _fname << ": n_variables = " << _n_variables << endl; + cout << "SPH_db_SpEC_H5_v2 " << _fname << ": lmax = " << _lmax << endl; + cout << "SPH_db_SpEC_H5_v2 " << _fname << ": timesteps = " << _iterations.size() << endl; + +} + + + + + +}
File [modified]: sph_database.hh Delta lines: +42 -1 =================================================================== --- PITTNullCode/SphericalHarmonicReconGen/src/sph_database.hh 2013-05-27 17:02:26 UTC (rev 84) +++ PITTNullCode/SphericalHarmonicReconGen/src/sph_database.hh 2013-06-04 20:40:27 UTC (rev 85) @@ -314,7 +314,7 @@
// the hdf5 file handle hid_t file; - + static const string time_attrib_name; static const string step_basename;
@@ -340,9 +340,50 @@ };
+/** + A SpEC-HDF5-v2 database +*/ +class SPH_db_SpEC_H5_v2 : public SPH_database +{ + public : + SPH_db_SpEC_H5_v2(const string& fname_, + const int cached_timesteps); + + virtual ~SPH_db_SpEC_H5_v2(); + + // scans the input file + virtual void scan(); + + // read all modes on all extraction spheres for a given timestep and variable number into array + virtual void read(const size_t timestep, const size_t varno, + vector<vector<vector<CCTK_COMPLEX> > >& coeff) const; + + private : + int remove_non_monotonic_steps(const int k); + void scan_HDF5(); + + + // table of varibale names: this maps a column number to a group. + // The first entry per varibale is the main group, the second is a possible subgroup, e.g. + // /g/Step00001/xx + // for SpEC + static const string varname_table[30]; + + string _fname; + + // the hdf5 file handle + hid_t file;
+ + + vector<CCTK_REAL> steps; + +};
+ + + /** This class represents a set of spherical coefficients for one variable and over all extraction spheres. It can read-in the coefficients from a database and reconstruct the angular dependence by
File [modified]: startup.cc Delta lines: +2 -0 =================================================================== --- PITTNullCode/SphericalHarmonicReconGen/src/startup.cc 2013-05-27 17:02:26 UTC (rev 84) +++ PITTNullCode/SphericalHarmonicReconGen/src/startup.cc 2013-06-04 20:40:27 UTC (rev 85) @@ -124,6 +124,8 @@ SHR::db[i] = new SHR::SPH_db_DAT(pathname+fnames[i], cached_timesteps+1, lmaxInFile, column_time, column_iteration, column_radius, column_lmax, column_n_variables, column_data); if (CCTK_EQUALS(format, "SpEC-H5")) SHR::db[i] = new SHR::SPH_db_SpEC_H5(pathname+fnames[i], cached_timesteps+1); + if (CCTK_EQUALS(format, "SpEC-H5-v2")) + SHR::db[i] = new SHR::SPH_db_SpEC_H5_v2(pathname+fnames[i], cached_timesteps+1);
// check if timestep-size found in file matches characteristic timestep! if (fabs(SHR::db[i]->delta_t()-CCTK_DELTA_TIME) > 1e-8)