Skip to content
Snippets Groups Projects
Commit 2389e7a0 authored by Giovanni La Mura's avatar Giovanni La Mura
Browse files

Implement HDF5 configuration binary output

parent 0bcc235f
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
src/cluster/Makefile
+ 7
− 4
View file @ 2389e7a0
......@@ -3,8 +3,8 @@ FC=gfortran
FCFLAGS=-std=legacy -O3
LFLAGS=
CXX=g++
CXXFLAGS=-O2 -ggdb -pg -coverage
CXXLFLAGS=
CXXFLAGS=-O3 -ggdb -pg -coverage
CXXLFLAGS=-L/usr/lib64 -lhdf5_hl -lhdf5
all: clu edfb np_cluster
......@@ -14,8 +14,8 @@ clu: clu.o
edfb: edfb.o
$(FC) $(FCFLAGS) -o $(BUILDDIR)/edfb $(BUILDDIR)/edfb.o $(LFLAGS)
np_cluster: $(BUILDDIR)/np_cluster.o $(BUILDDIR)/Commons.o $(BUILDDIR)/Configuration.o $(BUILDDIR)/Parsers.o $(BUILDDIR)/sph_subs.o $(BUILDDIR)/clu_subs.o $(BUILDDIR)/cluster.o
$(CXX) $(CXXFLAGS) $(CXXLFLAGS) -o $(BUILDDIR)/np_cluster $(BUILDDIR)/np_cluster.o $(BUILDDIR)/Commons.o $(BUILDDIR)/Configuration.o $(BUILDDIR)/Parsers.o $(BUILDDIR)/sph_subs.o $(BUILDDIR)/clu_subs.o $(BUILDDIR)/cluster.o
np_cluster: $(BUILDDIR)/np_cluster.o $(BUILDDIR)/Commons.o $(BUILDDIR)/Configuration.o $(BUILDDIR)/file_io.o $(BUILDDIR)/Parsers.o $(BUILDDIR)/sph_subs.o $(BUILDDIR)/clu_subs.o $(BUILDDIR)/cluster.o
$(CXX) $(CXXFLAGS) $(CXXLFLAGS) -o $(BUILDDIR)/np_cluster $(BUILDDIR)/np_cluster.o $(BUILDDIR)/Commons.o $(BUILDDIR)/Configuration.o $(BUILDDIR)/file_io.o $(BUILDDIR)/Parsers.o $(BUILDDIR)/sph_subs.o $(BUILDDIR)/clu_subs.o $(BUILDDIR)/cluster.o
$(BUILDDIR)/np_cluster.o:
$(CXX) $(CXXFLAGS) -c np_cluster.cpp -o $(BUILDDIR)/np_cluster.o
......@@ -26,6 +26,9 @@ $(BUILDDIR)/Commons.o:
$(BUILDDIR)/Configuration.o:
$(CXX) $(CXXFLAGS) -c ../libnptm/Configuration.cpp -o $(BUILDDIR)/Configuration.o
$(BUILDDIR)/file_io.o:
$(CXX) $(CXXFLAGS) -c ../libnptm/file_io.cpp -o $(BUILDDIR)/file_io.o
$(BUILDDIR)/Parsers.o:
$(CXX) $(CXXFLAGS) -c ../libnptm/Parsers.cpp -o $(BUILDDIR)/Parsers.o
......
src/cluster/cluster.cpp
+ 3
− 2
View file @ 2389e7a0
/*! \file cluster.cpp
*/
#include <cstdio>
#include <complex>
#include <exception>
#include <fstream>
#include <string>
#include <complex>
#ifndef INCLUDE_CONFIGURATION_H_
#include <exception>
#include "../include/Configuration.h"
#endif
......@@ -42,6 +42,7 @@ void cluster(string config_file, string data_file, string output_path) {
}
sconf->write_formatted(output_path + "/c_OEDFB");
sconf->write_binary(output_path + "/c_TEDF");
sconf->write_binary(output_path + "/c_TEDF.hd5", "HDF5");
GeometryConfiguration *gconf = NULL;
try {
gconf = GeometryConfiguration::from_legacy(data_file);
......
src/include/Configuration.h
+ 48
− 4
View file @ 2389e7a0
......@@ -218,6 +218,51 @@ protected:
//! \brief Flag to control whether to add an external layer.
bool use_external_sphere;
/*! \brief Build configuration from a HDF5 binary input file.
*
* This is the function called by the public method `from_binary()` in case of
* HDF5 mode selection. This function creates a configuration structure from
* a binary file written according to the HDF5 format standard.
*
* \param file_name: `string` Name of the binary configuration data file.
* \return config: `ScattererConfiguration*` Pointer to object containing the
* scatterer configuration data.
*/
static ScattererConfiguration *from_hdf5(std::string file_name);
/*! \brief Build configuration from legacy binary input file.
*
* This is the function called by the public method `from_binary()` in case of
* legacy mode selection. This function creates a configuration structure from
* a binary file written according to the proprietary mode used by the original
* FORTRAN code.
*
* \param file_name: `string` Name of the binary configuration data file.
* \return config: `ScattererConfiguration*` Pointer to object containing the
* scatterer configuration data.
*/
static ScattererConfiguration *from_legacy(std::string file_name);
/*! \brief Write the scatterer configuration data to HDF5 binary output.
*
* This function is invoked by the public method `write_binary()` with the
* "HDF5" format mode. It undertakes the task of writing the configuration
* information to a binary file using the standard HDF5 format.
*
* \param file_name: `string` Name of the binary configuration data file.
*/
void write_hdf5(std::string file_name);
/*! \brief Write the scatterer configuration data to legacy binary output.
*
* This function is invoked by the public method `write_binary()` with the
* "LEGACY" format mode. It undertakes the task of writing the configuration
* information to a binary file using a proprietary format, as it was done
* originally in the FORTRAN code.
*
* \param file_name: `string` Name of the binary configuration data file.
*/
void write_legacy(std::string file_name);
public:
/*! \brief Build a scatterer configuration structure.
*
......@@ -306,14 +351,13 @@ public:
* be useful to save the configuration data. `ScattererConfiguration.write_binary()`
* performs the operation of saving the configuration in binary format. This function
* can work in legacy mode, to write backward compatible configuration files, as well
* as by wrapping the data into common scientific formats (NB: this last option still
* needs to be implemented).
* as by wrapping the data into common scientific formats.
*
* \param file_name: `string` Name of the file to be written.
* \param mode: `string` Binary encoding. Can be one of "LEGACY", ... . Optional
* \param mode: `string` Binary encoding. Can be one of ["LEGACY", "HDF5"] . Optional
* (default is "LEGACY").
*/
void write_binary(std::string file_name, std::string mode = "LEGACY");
void write_binary(std::string file_name, std::string mode="LEGACY");
/*! \brief Write the scatterer configuration data to formatted text output.
*
......
src/include/file_io.h
+ 0
− 4
View file @ 2389e7a0
......@@ -120,10 +120,6 @@ class HDFFile {
*/
bool is_open() { return file_open_flag; }
/* ! \brief Read data from attached file.
*/
// herr_t read();
/*! \brief Write data to attached file.
*
* \param dataset_name: `string` Name of the dataset to write to.
......
src/libnptm/Configuration.cpp
+ 237
− 124
View file @ 2389e7a0
......@@ -4,9 +4,11 @@
#include <cmath>
#include <complex>
#include <cstdio>
#include <exception>
#include <fstream>
#include <regex>
#include <string>
#include "hdf5.h"
#ifndef INCLUDE_LIST_H_
#include "../include/List.h"
......@@ -20,6 +22,10 @@
#include "../include/Configuration.h"
#endif
#ifndef INCLUDE_FILE_IO_H_
#include "../include/file_io.h"
#endif
using namespace std;
GeometryConfiguration::GeometryConfiguration(
......@@ -203,8 +209,8 @@ ScattererConfiguration::ScattererConfiguration(
radii_of_spheres = ros_vector;
nshl_vec = nshl_vector;
rcf = rcf_vector;
idfc = dielectric_func_type,
dc0_matrix = dc_matrix;
idfc = dielectric_func_type;
dc0_matrix = dc_matrix;
use_external_sphere = is_external;
exdc = ex;
wp = w;
......@@ -248,6 +254,24 @@ ScattererConfiguration::~ScattererConfiguration() {
}
ScattererConfiguration* ScattererConfiguration::from_binary(string file_name, string mode) {
ScattererConfiguration *conf = NULL;
if (mode.compare("LEGACY") == 0) {
conf = ScattererConfiguration::from_legacy(file_name);
} else if (mode.compare("HDF5") == 0) {
conf = ScattererConfiguration::from_hdf5(file_name);
} else {
string message = "Unknown format mode: \"" + mode + "\"";
throw UnrecognizedConfigurationException(message);
}
return conf;
}
ScattererConfiguration* ScattererConfiguration::from_hdf5(string file_name) {
ScattererConfiguration *conf = NULL;
return conf;
}
ScattererConfiguration* ScattererConfiguration::from_legacy(string file_name) {
int nsph;
int *iog;
double _exdc, _wp, _xip;
......@@ -257,80 +281,78 @@ ScattererConfiguration* ScattererConfiguration::from_binary(string file_name, st
double *ros_vector;
double **rcf_vector;
complex<double> ***dc0m;
if (mode.compare("LEGACY") == 0) { // Legacy mode was chosen.
int max_ici = 0;
fstream input;
input.open(file_name.c_str(), ios::in | ios::binary);
if (input.is_open()) {
input.read(reinterpret_cast<char *>(&nsph), sizeof(int));
iog = new int[nsph]();
for (int i = 0; i < nsph; i++) {
input.read(reinterpret_cast<char *>(&(iog[i])), sizeof(int));
int max_ici = 0;
fstream input;
input.open(file_name.c_str(), ios::in | ios::binary);
if (input.is_open()) {
input.read(reinterpret_cast<char *>(&nsph), sizeof(int));
iog = new int[nsph]();
for (int i = 0; i < nsph; i++) {
input.read(reinterpret_cast<char *>(&(iog[i])), sizeof(int));
}
input.read(reinterpret_cast<char *>(&_exdc), sizeof(double));
input.read(reinterpret_cast<char *>(&_wp), sizeof(double));
input.read(reinterpret_cast<char *>(&_xip), sizeof(double));
input.read(reinterpret_cast<char *>(&_idfc), sizeof(int));
input.read(reinterpret_cast<char *>(&nxi), sizeof(int));
try {
xi_vec = new double[nxi]();
} catch (const bad_alloc &ex) {
throw UnrecognizedConfigurationException("Wrong parameter set: invalid number of scales " + to_string(nxi));
}
for (int i = 0; i < nxi; i++) {
input.read(reinterpret_cast<char *>(&(xi_vec[i])), sizeof(double));
}
nshl_vector = new int[nsph]();
ros_vector = new double[nsph]();
rcf_vector = new double*[nsph];
for (int i115 = 1; i115 <= nsph; i115++) {
if (iog[i115 - 1] < i115) {
rcf_vector[i115 - 1] = new double[1]();
continue;
}
input.read(reinterpret_cast<char *>(&_exdc), sizeof(double));
input.read(reinterpret_cast<char *>(&_wp), sizeof(double));
input.read(reinterpret_cast<char *>(&_xip), sizeof(double));
input.read(reinterpret_cast<char *>(&_idfc), sizeof(int));
input.read(reinterpret_cast<char *>(&nxi), sizeof(int));
input.read(reinterpret_cast<char *>(&(nshl_vector[i115 - 1])), sizeof(int));
input.read(reinterpret_cast<char *>(&(ros_vector[i115 - 1])), sizeof(double));
int nsh = nshl_vector[i115 - 1];
if (max_ici < (nsh + 1) / 2) max_ici = (nsh + 1) / 2;
try {
xi_vec = new double[nxi]();
rcf_vector[i115 - 1] = new double[nsh]();
} catch (const bad_alloc &ex) {
throw UnrecognizedConfigurationException("Wrong parameter set: invalid number of scales " + to_string(nxi));
throw UnrecognizedConfigurationException("Wrong parameter set: invalid number of layers " + to_string(nsh));
}
for (int i = 0; i < nxi; i++) {
input.read(reinterpret_cast<char *>(&(xi_vec[i])), sizeof(double));
for (int nsi = 0; nsi < nsh; nsi++) {
input.read(reinterpret_cast<char *>(&(rcf_vector[i115 - 1][nsi])), sizeof(double));
}
nshl_vector = new int[nsph]();
ros_vector = new double[nsph]();
rcf_vector = new double*[nsph];
for (int i115 = 1; i115 <= nsph; i115++) {
if (iog[i115 - 1] < i115) {
rcf_vector[i115 - 1] = new double[1]();
continue;
}
input.read(reinterpret_cast<char *>(&(nshl_vector[i115 - 1])), sizeof(int));
input.read(reinterpret_cast<char *>(&(ros_vector[i115 - 1])), sizeof(double));
int nsh = nshl_vector[i115 - 1];
if (max_ici < (nsh + 1) / 2) max_ici = (nsh + 1) / 2;
try {
rcf_vector[i115 - 1] = new double[nsh]();
} catch (const bad_alloc &ex) {
throw UnrecognizedConfigurationException("Wrong parameter set: invalid number of layers " + to_string(nsh));
}
for (int nsi = 0; nsi < nsh; nsi++) {
input.read(reinterpret_cast<char *>(&(rcf_vector[i115 - 1][nsi])), sizeof(double));
}
}
dc0m = new complex<double>**[max_ici];
for (int dim1 = 0; dim1 < max_ici; dim1++) {
dc0m[dim1] = new complex<double>*[nsph];
for (int dim2 = 0; dim2 < nsph; dim2++) {
dc0m[dim1][dim2] = new complex<double>[nxi]();
}
}
dc0m = new complex<double>**[max_ici];
int dim3 = (_idfc == 0) ? nxi : 1;
for (int dim1 = 0; dim1 < max_ici; dim1++) {
dc0m[dim1] = new complex<double>*[nsph];
for (int dim2 = 0; dim2 < nsph; dim2++) {
dc0m[dim1][dim2] = new complex<double>[dim3]();
}
for (int jxi468 = 1; jxi468 <= nxi; jxi468++) {
if (_idfc != 0 && jxi468 > 1) continue;
for (int i162 = 1; i162 <= nsph; i162++) {
if (iog[i162 - 1] < i162) continue;
int nsh = nshl_vector[i162 - 1];
int ici = (nsh + 1) / 2; // QUESTION: is integer division really intended here?
for (int i157 = 0; i157 < ici; i157++) {
double dc0_real, dc0_img;
input.read(reinterpret_cast<char *>(&dc0_real), sizeof(double));
input.read(reinterpret_cast<char *>(&dc0_img), sizeof(double));
dc0m[i157][i162 - 1][jxi468 - 1] = dc0_real + 1i * dc0_img;
}
}
for (int jxi468 = 1; jxi468 <= nxi; jxi468++) {
if (_idfc != 0 && jxi468 > 1) continue;
for (int i162 = 1; i162 <= nsph; i162++) {
if (iog[i162 - 1] < i162) continue;
int nsh = nshl_vector[i162 - 1];
int ici = (nsh + 1) / 2; // QUESTION: is integer division really intended here?
for (int i157 = 0; i157 < ici; i157++) {
double dc0_real, dc0_img;
input.read(reinterpret_cast<char *>(&dc0_real), sizeof(double));
input.read(reinterpret_cast<char *>(&dc0_img), sizeof(double));
dc0m[i157][i162 - 1][jxi468 - 1] = dc0_real + 1i * dc0_img;
}
}
input.close();
} else { // Opening of the input file did not succeed.
OpenConfigurationFileException ex(file_name);
throw ex;
}
} else { // A different binary format was chosen.
//TODO: this part is not yet implemented.
// Functions to write optimized file formats may be invoked here.
input.close();
} else { // Opening of the input file did not succeed.
OpenConfigurationFileException ex(file_name);
throw ex;
}
ScattererConfiguration *conf = new ScattererConfiguration(
nsph,
xi_vec,
......@@ -511,20 +533,10 @@ ScattererConfiguration* ScattererConfiguration::from_dedfb(string dedfb_file_nam
}
}
}
//last_read_line++;
int *iog_vector = new int[nsph]();
double *ros_vector = new double[nsph]();
double **rcf_vector = new double*[nsph];
int *nshl_vector = new int[nsph]();
/*for (int i = 0; i < nsph; i++) {
string read_format = "";
for (int j = 0; j < (i % 15); j++) read_format += " %*d";
read_format += " %d";
sscanf(file_lines[last_read_line].c_str(), read_format.c_str(), (iog_vector + i));
if (i > 0 && i % 15 == 0) {
last_read_line++;
}
}*/
int filled_iogs = 0;
re = regex("[0-9]+");
while (filled_iogs < nsph) {
......@@ -561,8 +573,6 @@ ScattererConfiguration* ScattererConfiguration::from_dedfb(string dedfb_file_nam
rcf_vector[i113 - 1] = new double[nsh]();
for (int ns = 0; ns < nsh; ns++) {
double ns_rcf;
//int ns_rcf_exp;
//sscanf(file_lines[++last_read_line].c_str(), " %lf D%d", &ns_rcf, &ns_rcf_exp);
str_target = file_lines[++last_read_line];
regex_search(str_target, m, re);
str_number = m.str();
......@@ -573,10 +583,11 @@ ScattererConfiguration* ScattererConfiguration::from_dedfb(string dedfb_file_nam
}
}
complex<double> ***dc0m = new complex<double>**[max_ici];
int dim3 = (_idfc == 0) ? nxi : 1;
for (int dim1 = 0; dim1 < max_ici; dim1++) {
dc0m[dim1] = new complex<double>*[nsph];
for (int dim2 = 0; dim2 < nsph; dim2++) {
dc0m[dim1][dim2] = new complex<double>[nxi]();
dc0m[dim1][dim2] = new complex<double>[dim3]();
}
}
re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?[0-9]+)?");
......@@ -589,8 +600,6 @@ ScattererConfiguration* ScattererConfiguration::from_dedfb(string dedfb_file_nam
if (i162 == 1) ici = ici + ies;
for (int i157 = 0; i157 < ici; i157++) {
double dc0_real, dc0_img;
//int dc0_real_exp, dc0_img_exp;
//sscanf(file_lines[++last_read_line].c_str(), " (%lf D%d, %lf D%d)", &dc0_real, &dc0_real_exp, &dc0_img, &dc0_img_exp);
str_target = file_lines[++last_read_line];
regex_search(str_target, m, re);
string str_number = m.str();
......@@ -603,7 +612,7 @@ ScattererConfiguration* ScattererConfiguration::from_dedfb(string dedfb_file_nam
str_number = regex_replace(str_number, regex("D"), "e");
str_number = regex_replace(str_number, regex("d"), "e");
dc0_img = stod(str_number);
dc0m[i157][i162 - 1][jxi468 - 1] = std::complex<double>(dc0_real, dc0_img);
dc0m[i157][i162 - 1][jxi468 - 1] = complex<double>(dc0_real, dc0_img);
}
}
}
......@@ -676,55 +685,159 @@ void ScattererConfiguration::print() {
}
void ScattererConfiguration::write_binary(string file_name, string mode) {
if (mode.compare("LEGACY") == 0) {
write_legacy(file_name);
} else if (mode.compare("HDF5") == 0) {
write_hdf5(file_name);
} else {
string message = "Unknown format mode: \"" + mode + "\"";
throw UnrecognizedConfigurationException(message);
}
}
void ScattererConfiguration::write_hdf5(string file_name) {
const double two_pi = acos(0.0) * 4.0;
const double evc = 6.5821188e-16;
int ies = (use_external_sphere)? 1 : 0;
int max_ici = 0;
bool is_new_vector = false;
if (mode.compare("LEGACY") == 0) { // Legacy mode was chosen.
fstream output;
int ies = (use_external_sphere)? 1 : 0;
output.open(file_name.c_str(), ios::out | ios::binary);
output.write(reinterpret_cast<char *>(&number_of_spheres), sizeof(int));
for (int i = 0; i < number_of_spheres; i++)
output.write(reinterpret_cast<char *>(&(iog_vec[i])), sizeof(int));
output.write(reinterpret_cast<char *>(&exdc), sizeof(double));
output.write(reinterpret_cast<char *>(&wp), sizeof(double));
output.write(reinterpret_cast<char *>(&xip), sizeof(double));
output.write(reinterpret_cast<char *>(&idfc), sizeof(int));
output.write(reinterpret_cast<char *>(&number_of_scales), sizeof(int));
for (int i = 0; i < number_of_scales; i++)
output.write(reinterpret_cast<char *>(&(scale_vec[i])), sizeof(double));
for (int i115 = 1; i115 <= number_of_spheres; i115++) {
if (iog_vec[i115 - 1] < i115) continue;
output.write(reinterpret_cast<char *>(&(nshl_vec[i115 - 1])), sizeof(int));
output.write(reinterpret_cast<char *>(&(radii_of_spheres[i115 - 1])), sizeof(double));
int nsh = nshl_vec[i115 - 1];
if (i115 == 1) nsh += ies;
if (max_ici < (nsh + 1) / 2) max_ici = (nsh + 1) / 2;
for (int nsi = 0; nsi < nsh; nsi++)
output.write(reinterpret_cast<char *>(&(rcf[i115 - 1][nsi])), sizeof(double));
List<string> rec_name_list(1);
List<string> rec_type_list(1);
List<void *> rec_ptr_list(1);
string str_type, str_name;
rec_name_list.set(0, "NSPH");
rec_type_list.set(0, "INT32_(1)");
rec_ptr_list.set(0, &number_of_spheres);
rec_name_list.append("IOGVEC");
str_type = "INT32_(" + to_string(number_of_spheres) + ")";
rec_type_list.append(str_type);
rec_ptr_list.append(iog_vec);
rec_name_list.append("EXDC");
rec_type_list.append("FLOAT64_(1)");
rec_ptr_list.append(&exdc);
rec_name_list.append("WP");
rec_type_list.append("FLOAT64_(1)");
rec_ptr_list.append(&wp);
rec_name_list.append("XIP");
rec_type_list.append("FLOAT64_(1)");
rec_ptr_list.append(&xip);
rec_name_list.append("IDFC");
rec_type_list.append("INT32_(1)");
rec_ptr_list.append(&idfc);
rec_name_list.append("NXI");
rec_type_list.append("INT32_(1)");
rec_ptr_list.append(&number_of_scales);
rec_name_list.append("XIVEC");
str_type = "FLOAT64_(" + to_string(number_of_scales) + ")";
rec_type_list.append(str_type);
rec_ptr_list.append(scale_vec);
for (int i115 = 1; i115 <= number_of_spheres; i115++) {
if (iog_vec[i115 - 1] < i115) continue;
str_name = "NSHL_" + to_string(i115);
rec_name_list.append(str_name);
rec_type_list.append("INT32_(1)");
rec_ptr_list.append(&(nshl_vec[i115 - 1]));
str_name = "ROS_" + to_string(i115);
rec_name_list.append(str_name);
rec_type_list.append("FLOAT64_(1)");
rec_ptr_list.append(&(radii_of_spheres[i115 - 1]));
int nsh = nshl_vec[i115 - 1];
if (i115 == 1) nsh += ies;
if (max_ici < (nsh + 1) / 2) max_ici = nsh + 1 / 2;
str_name = "RCF_" + to_string(i115);
str_type = "FLOAT64_(" + to_string(nsh) + ")";
rec_name_list.append(str_name);
rec_type_list.append(str_type);
rec_ptr_list.append(&(rcf[i115 - 1][0]));
}
int dim3 = (idfc == 0) ? number_of_scales : 1;
int dc0m_size = 2 * dim3 * number_of_spheres * max_ici;
double *dc0m = new double[dc0m_size];
int dc0_index = 0;
for (int jxi468 = 1; jxi468 <= number_of_scales; jxi468++) {
if (idfc != 0 && jxi468 > 1) continue;
for (int i162 = 1; i162 <= number_of_spheres; i162++) {
if (iog_vec[i162 - 1] < i162) continue;
int nsh = nshl_vec[i162 - 1];
int ici = (nsh + 1) / 2;
if (i162 == 1) ici = ici + ies;
for (int i157 = 0; i157 < ici; i157++) {
double dc0_real, dc0_imag;
dc0_real = dc0_matrix[i157][i162 - 1][jxi468 - 1].real();
dc0_imag = dc0_matrix[i157][i162 - 1][jxi468 - 1].imag();
dc0m[dc0_index++] = dc0_real;
dc0m[dc0_index++] = dc0_imag;
}
}
for (int jxi468 = 1; jxi468 <= number_of_scales; jxi468++) {
if (idfc != 0 && jxi468 > 1) continue;
for (int i162 = 1; i162 <= number_of_spheres; i162++) {
if (iog_vec[i162 - 1] < i162) continue;
int nsh = nshl_vec[i162 - 1];
int ici = (nsh + 1) / 2; // QUESTION: is integer division really intended here?
if (i162 == 1) ici = ici + ies;
for (int i157 = 0; i157 < ici; i157++) {
double dc0_real, dc0_img;
dc0_real = dc0_matrix[i157][i162 - 1][jxi468 - 1].real();
dc0_img = dc0_matrix[i157][i162 - 1][jxi468 - 1].imag();
// The FORTRAN code writes the complex numbers as a 16-byte long binary stream.
// Here we assume that the 16 bytes are equally split in 8 bytes to represent the
// real part and 8 bytes to represent the imaginary one.
output.write(reinterpret_cast<char *>(&dc0_real), sizeof(double));
output.write(reinterpret_cast<char *>(&dc0_img), sizeof(double));
}
}
str_type = "FLOAT64_(" + to_string(dc0m_size) + ")";
rec_name_list.append("DC0M");
rec_type_list.append(str_type);
rec_ptr_list.append(dc0m);
string *rec_names = rec_name_list.to_array();
string *rec_types = rec_type_list.to_array();
void **rec_pointers = rec_ptr_list.to_array();
const int rec_num = rec_name_list.length();
FileSchema schema(rec_num, rec_types, rec_names);
HDFFile *hdf_file = HDFFile::from_schema(schema, "c_TEDF.hd5", H5F_ACC_TRUNC);
for (int ri = 0; ri < rec_num; ri++)
hdf_file->write(rec_names[ri], rec_types[ri], rec_pointers[ri]);
hdf_file->close();
// Clean memory
delete[] dc0m;
delete[] rec_names;
delete[] rec_types;
delete[] rec_pointers;
delete hdf_file;
}
void ScattererConfiguration::write_legacy(string file_name) {
const double two_pi = acos(0.0) * 4.0;
const double evc = 6.5821188e-16;
fstream output;
int ies = (use_external_sphere)? 1 : 0;
output.open(file_name.c_str(), ios::out | ios::binary);
output.write(reinterpret_cast<char *>(&number_of_spheres), sizeof(int));
for (int i = 0; i < number_of_spheres; i++)
output.write(reinterpret_cast<char *>(&(iog_vec[i])), sizeof(int));
output.write(reinterpret_cast<char *>(&exdc), sizeof(double));
output.write(reinterpret_cast<char *>(&wp), sizeof(double));
output.write(reinterpret_cast<char *>(&xip), sizeof(double));
output.write(reinterpret_cast<char *>(&idfc), sizeof(int));
output.write(reinterpret_cast<char *>(&number_of_scales), sizeof(int));
for (int i = 0; i < number_of_scales; i++)
output.write(reinterpret_cast<char *>(&(scale_vec[i])), sizeof(double));
for (int i115 = 1; i115 <= number_of_spheres; i115++) {
if (iog_vec[i115 - 1] < i115) continue;
output.write(reinterpret_cast<char *>(&(nshl_vec[i115 - 1])), sizeof(int));
output.write(reinterpret_cast<char *>(&(radii_of_spheres[i115 - 1])), sizeof(double));
int nsh = nshl_vec[i115 - 1];
if (i115 == 1) nsh += ies;
for (int nsi = 0; nsi < nsh; nsi++)
output.write(reinterpret_cast<char *>(&(rcf[i115 - 1][nsi])), sizeof(double));
}
for (int jxi468 = 1; jxi468 <= number_of_scales; jxi468++) {
if (idfc != 0 && jxi468 > 1) continue;
for (int i162 = 1; i162 <= number_of_spheres; i162++) {
if (iog_vec[i162 - 1] < i162) continue;
int nsh = nshl_vec[i162 - 1];
int ici = (nsh + 1) / 2; // QUESTION: is integer division really intended here?
if (i162 == 1) ici = ici + ies;
for (int i157 = 0; i157 < ici; i157++) {
double dc0_real, dc0_img;
dc0_real = dc0_matrix[i157][i162 - 1][jxi468 - 1].real();
dc0_img = dc0_matrix[i157][i162 - 1][jxi468 - 1].imag();
// The FORTRAN code writes the complex numbers as a 16-byte long binary stream.
// Here we assume that the 16 bytes are equally split in 8 bytes to represent the
// real part and 8 bytes to represent the imaginary one.
output.write(reinterpret_cast<char *>(&dc0_real), sizeof(double));
output.write(reinterpret_cast<char *>(&dc0_img), sizeof(double));
}
}
output.close();
}
output.close();
}
void ScattererConfiguration::write_formatted(string file_name) {
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment