User: reisswig Date: 2011/11/07 06:03 PM

Modified: /PITTNullCode/SphericalHarmonicReconASCII/ configuration.ccl, param.ccl /PITTNullCode/SphericalHarmonicReconASCII/src/ sph_database.cc, sph_database.hh, startup.cc

Log: Add capability of hdf5.

File Changes:

Directory: /PITTNullCode/SphericalHarmonicReconASCII/ =====================================================

File [modified]: configuration.ccl Delta lines: +1 -1 =================================================================== --- PITTNullCode/SphericalHarmonicReconASCII/configuration.ccl 2011-09-04 20:16:09 UTC (rev 56) +++ PITTNullCode/SphericalHarmonicReconASCII/configuration.ccl 2011-11-08 00:03:19 UTC (rev 57) @@ -1 +1 @@ -REQUIRES +REQUIRES HDF5

File [modified]: param.ccl Delta lines: +3 -2 =================================================================== --- PITTNullCode/SphericalHarmonicReconASCII/param.ccl 2011-09-04 20:16:09 UTC (rev 56) +++ PITTNullCode/SphericalHarmonicReconASCII/param.ccl 2011-11-08 00:03:19 UTC (rev 57) @@ -5,8 +5,9 @@

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 :: "all variables and lm-modes are stored in one row for a given timestep" + 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"

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/SphericalHarmonicReconASCII/src/ =========================================================

File [modified]: sph_database.cc Delta lines: +227 -0 =================================================================== --- PITTNullCode/SphericalHarmonicReconASCII/src/sph_database.cc 2011-09-04 20:16:09 UTC (rev 56) +++ PITTNullCode/SphericalHarmonicReconASCII/src/sph_database.cc 2011-11-08 00:03:19 UTC (rev 57) @@ -855,7 +855,234 @@ }

+ + + + + + + + + + + + + +//------------------------------------------------------------------------------------------- +// +// +// SpEC-H5 format goes here........... +// +// +//------------------------------------------------------------------------------------------ + + + + +const string SPH_db_SpEC_H5::time_attrib_name = "Time"; +const string SPH_db_SpEC_H5::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::varname_table[30][2] = + { { "Lapse", "scalar"}, + { "Shift", "x" }, { "Shift", "y" }, { "Shift", "z" }, + { "g", "xx" }, { "g", "xy" }, { "g", "xz" }, + { "g", "yy" }, { "g", "yz" }, + { "g", "zz" }, + + { "DrLapse", "scalar"}, + { "DrShift", "x" }, { "DrShift", "y" }, { "DrShift", "z" }, + { "Drg", "xx" }, { "Drg", "xy" }, { "Drg", "xz" }, + { "Drg", "yy" }, { "Drg", "yz" }, + { "Drg", "zz" }, + + { "DtLapse", "scalar"}, + { "DtShift", "x" }, { "DtShift", "y" }, { "DtShift", "z" }, + { "Dtg", "xx" }, { "Dtg", "xy" }, { "Dtg", "xz" }, + { "Dtg", "yy" }, { "Dtg", "yz" }, + { "Dtg", "zz" } }; + + +SPH_db_SpEC_H5::SPH_db_SpEC_H5(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::~SPH_db_SpEC_H5() +{ + // close file + H5Fclose(file); +} +

+int SPH_db_SpEC_H5::remove_non_monotonic_steps(const int k) +{ + // go backwards from k and remove until we get monotonicity + const double time = steps[k].time; + const int size = steps.size(); + int c = 1; + for (int i=k-1; i >= 0; --i) { + if (time > steps[i].time) { + 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::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 (int 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::read(const int timestep, const int 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][0]+"/"+steps[timestep].name+"/"+varname_table[varno][1]).c_str());//, H5P_DEFAULT); + vector<double> buffer((_lmax+1)*(_lmax+1)*2); + H5Dread(dataset, H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, H5P_DEFAULT, &buffer.front()); + 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]); + } + } +} + + + +herr_t SPH_db_SpEC_H5::H5iter(hid_t loc_id, const char* name, const H5L_info_t* info, void* operator_data) +{ + SPH_db_SpEC_H5* SpEC_H5 = (SPH_db_SpEC_H5*) operator_data; + + // Datasets are ignored! We only iterate over groups + switch (info->type) { + case H5O_TYPE_GROUP: { + const string gname = name; + // we have found a timstep group if the group contains the "step_basename" substring + if (gname.find(step_basename) != string::npos) { + // get time attribute + double time; + hid_t group = H5Gopen(loc_id, name); + hid_t attribute = H5Aopen(group, time_attrib_name.c_str(), H5P_DEFAULT); + H5Aread(attribute, H5T_NATIVE_DOUBLE, &time); + H5Aclose(attribute); + H5Gclose(group); + + SpEC_H5->steps.push_back(timestep_t(name, time)); + } + break; + } + case H5O_TYPE_DATASET: + // ignore datasets + return 0; + break; + default: + break; + } + + return 0; +} + + + +void SPH_db_SpEC_H5::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) + herr_t status = 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 (int 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, (varname_table[0][0]+"/"+steps[0].name+"/"+varname_table[0][1]).c_str()); //, H5P_DEFAULT); + //cout << (varname_table[0][0]+"/"+steps[0].name+"/"+varname_table[0][1]) << endl; + hid_t dataspace = H5Dget_space(dataset); + hsize_t size = H5Sget_simple_extent_npoints(dataspace); + H5Sclose(dataspace); + H5Dclose(dataset); + + _lmax = sqrt(size/2)-1; + _n_variables = 30; + _n_spheres = 1; + + // cross-check that all expected variables are present + + + // display some info + cout << "SPH_db_SpEC_H5 " << _fname << ": n_spheres = " << _n_spheres << endl; + cout << "SPH_db_SpEC_H5 " << _fname << ": n_variables = " << _n_variables << endl; + cout << "SPH_db_SpEC_H5 " << _fname << ": lmax = " << _lmax << endl; + cout << "SPH_db_SpEC_H5 " << _fname << ": timesteps = " << _iterations.size() << endl; +} + + + + + + + +}

File [modified]: sph_database.hh Delta lines: +63 -0 =================================================================== --- PITTNullCode/SphericalHarmonicReconASCII/src/sph_database.hh 2011-09-04 20:16:09 UTC (rev 56) +++ PITTNullCode/SphericalHarmonicReconASCII/src/sph_database.hh 2011-11-08 00:03:19 UTC (rev 57) @@ -18,10 +18,14 @@

#include "mpi.h"

+#define H5_USE_16_API 1 +#include "hdf5.h"

#ifndef _SPH_database_ #define _SPH_database_

+ + namespace SHR {

using namespace std; @@ -276,7 +280,66 @@

+ /** + A SpEC-HDF5 database +*/ +class SPH_db_SpEC_H5 : public SPH_database +{ + public : + SPH_db_SpEC_H5(const string& fname_, + const int cached_timesteps); + + virtual ~SPH_db_SpEC_H5(); + + // 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 int timestep, const int varno, + vector<vector<vector<CCTK_COMPLEX> > >& coeff) const; + + private : + + int remove_non_monotonic_steps(const int k); + void scan_HDF5(); + static herr_t H5iter(hid_t loc_id, const char* name, const H5L_info_t* info, void* operator_data); + + + string _fname; + + // the hdf5 file handle + hid_t file; + + static const string time_attrib_name; + static const string step_basename; + + // 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][2]; + + struct timestep_t + { + string name; + double time; + timestep_t(const string& name_, const double time_) + : name(name_), time(time_) + {} + bool operator<(const timestep_t& rhs) const { return time < rhs.time; } + }; + + /// the cylces and names of the steps contained in the file + vector<timestep_t> 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 recomposing with the basis-functions, i.e. the sYlm.

File [modified]: startup.cc Delta lines: +2 -0 =================================================================== --- PITTNullCode/SphericalHarmonicReconASCII/src/startup.cc 2011-09-04 20:16:09 UTC (rev 56) +++ PITTNullCode/SphericalHarmonicReconASCII/src/startup.cc 2011-11-08 00:03:19 UTC (rev 57) @@ -122,6 +122,8 @@ SHR::db[i] = new SHR::SPH_db_ASCII(pathname+fnames[i], cached_timesteps+1, column_time, column_iteration, column_radius); if (CCTK_EQUALS(format, "DAT")) 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);

// check if timestep-size found in file matches characteristic timestep! if (fabs(SHR::db[i]->delta_t()-CCTK_DELTA_TIME) > 1e-8)