diff --git a/src/libnptm/Configuration.cpp b/src/libnptm/Configuration.cpp
index 63faf394b7cc76573a1c248e0fde060d2b71a619..3485d543bef862afd5378ee32f6d061a9ba291ea 100644
--- a/src/libnptm/Configuration.cpp
+++ b/src/libnptm/Configuration.cpp
@@ -234,21 +234,21 @@ ScattererConfiguration::ScattererConfiguration(
}
ScattererConfiguration::~ScattererConfiguration() {
- int max_ici = 0;
+ int configurations = 0;
for (int i = 1; i <= number_of_spheres; i++) {
if (iog_vec[i - 1] < i) continue;
- int nsh = nshl_vec[i - 1];
- if (max_ici < (nsh + 1) / 2) max_ici = (nsh + 1) / 2;
+ configurations++;
}
- for (int i = 0; i < max_ici; i++) {
+ for (int i = 0; i < configurations; i++) {
for (int j = 0; j < number_of_spheres; j++) {
delete[] dc0_matrix[i][j];
}
}
delete[] dc0_matrix;
- for (int i = 0; i < number_of_spheres; i++) {
+ for (int i = 0; i < configurations; i++) {
delete[] rcf[i];
}
+ delete[] rcf;
delete[] nshl_vec;
delete[] radii_of_spheres;
delete[] iog_vec;
@@ -430,88 +430,82 @@ ScattererConfiguration* ScattererConfiguration::from_dedfb(string dedfb_file_nam
}
}
int *iog_vector = new int[nsph]();
- double *ros_vector = new double[nsph]();
- double **rcf_vector = new double*[nsph];
- int *nshl_vector = new int[nsph]();
- int filled_iogs = 0;
+ str_target = file_lines[++last_read_line];
re = regex("[0-9]+");
- while (filled_iogs < nsph) {
- string str_target = file_lines[++last_read_line];
- while(regex_search(str_target, m, re)) {
- iog_vector[filled_iogs++] = stoi(m.str());
+ int configurations = 0;
+ for (int i = 1; i <= nsph; i++) {
+ bool success = regex_search(str_target, m, re);
+ if (success) {
+ iog_vector[i - 1] = stoi(m.str());
str_target = m.suffix().str();
- if (filled_iogs == nsph) break;
+ if (iog_vector[i - 1] >= i) configurations++;
+ } else {
+ str_target = file_lines[++last_read_line];
+ i--;
}
}
+ int index = 0;
+ double *ros_vector = new double[configurations]();
+ int *nshl_vector = new int[configurations]();
+ double **rcf_vector = new double*[configurations];
for (int i113 = 1; i113 <= nsph; i113++) {
- int i_val, nsh;
- double ros_val;
- if (iog_vector[i113 - 1] < i113) {
- rcf_vector[i113 - 1] = new double[1]();
- continue;
- }
- re = regex("[0-9]+");
+ if (iog_vector[i113 - 1] < i113) continue;
str_target = file_lines[++last_read_line];
regex_search(str_target, m, re);
- i_val = stoi(m.str());
- str_target = m.suffix();
+ nshl_vector[i113 - 1] = stoi(m.str());
+ str_target = m.suffix().str();
re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?[0-9]+)?");
regex_search(str_target, m, re);
string str_number = m.str();
str_number = regex_replace(str_number, regex("D"), "e");
str_number = regex_replace(str_number, regex("d"), "e");
- ros_val = stod(str_number);
- nshl_vector[i113 - 1] = i_val;
- if (max_ici < (i_val + 1) / 2) max_ici = (i_val + 1) / 2;
- ros_vector[i113 - 1] = ros_val;
- nsh = nshl_vector[i113 - 1];
+ ros_vector[i113 - 1] = stod(m.str());
+ int nsh = nshl_vector[i113 - 1];
if (i113 == 1) nsh += ies;
rcf_vector[i113 - 1] = new double[nsh]();
- for (int ns = 0; ns < nsh; ns++) {
- double ns_rcf;
+ for (int ns112 = 0; ns112 < nsh; ns112++) {
str_target = file_lines[++last_read_line];
regex_search(str_target, m, re);
str_number = m.str();
str_number = regex_replace(str_number, regex("D"), "e");
str_number = regex_replace(str_number, regex("d"), "e");
- ns_rcf = stod(str_number);
- rcf_vector[i113 - 1][ns] = ns_rcf;
- }
- }
- 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>[dim3]();
- }
- }
- re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?[0-9]+)?");
+ rcf_vector[i113 - 1][ns112] = stod(str_number);
+ } // ns112 loop
+ } // i113 loop
+
+ complex<double> ***dc0m = new complex<double>**[configurations];
+ complex<double> *dc0 = new complex<double>[configurations]();
+ for (int di = 0; di < configurations; di++) {
+ dc0m[di] = new complex<double>*[nsph];
+ for (int dj = 0; dj < nsph; dj++) dc0m[di][dj] = new complex<double>[nxi]();
+ } // di loop
for (int jxi468 = 1; jxi468 <= nxi; jxi468++) {
- if (_idfc != 0 && jxi468 > 1) continue;
+ if (_idfc != 0 && jxi468 > 1) continue; // jxi468 loop
for (int i162 = 1; i162 <= nsph; i162++) {
- if (iog_vector[i162 - 1] < i162) continue;
- int nsh = nshl_vector[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;
- str_target = file_lines[++last_read_line];
- regex_search(str_target, m, re);
- string str_number = m.str();
- str_number = regex_replace(str_number, regex("D"), "e");
- str_number = regex_replace(str_number, regex("d"), "e");
- dc0_real = stod(str_number);
- str_target = m.suffix().str();
- regex_search(str_target, m, re);
- str_number = m.str();
- 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] = complex<double>(dc0_real, dc0_img);
+ if (iog_vector[i162 - 1] >= i162) {
+ int nsh = nshl_vector[i162 - 1];
+ int ici = (nsh + 1) / 2;
+ if (i162 == 1) ici += ies;
+ for (int ic157 = 0; ic157 < ici; ic157++) {
+ str_target = file_lines[++last_read_line];
+ regex_search(str_target, m, re);
+ string str_number = m.str();
+ str_number = regex_replace(str_number, regex("D"), "e");
+ str_number = regex_replace(str_number, regex("d"), "e");
+ double rval = stod(str_number);
+ str_target = m.suffix().str();
+ regex_search(str_target, m, re);
+ str_number = m.str();
+ str_number = regex_replace(str_number, regex("D"), "e");
+ str_number = regex_replace(str_number, regex("d"), "e");
+ double ival = stod(str_number);
+ dc0[ic157] = complex<double>(rval, ival);
+ dc0m[ic157][i162 - 1][jxi468 - 1] = dc0[ic157];
+ } // ic157 loop
}
- }
- }
+ } // i162 loop
+ } // jxi468 loop
+ delete[] dc0;
ScattererConfiguration *config = new ScattererConfiguration(
nsph,
@@ -627,104 +621,83 @@ ScattererConfiguration* ScattererConfiguration::from_hdf5(string file_name) {
}
ScattererConfiguration* ScattererConfiguration::from_legacy(string file_name) {
- int nsph;
- int *iog;
+ int _nsph;
+ int *_iog_vec;
+ int _ies;
double _exdc, _wp, _xip;
- int _idfc, nxi;
- int *nshl_vector;
- double *xi_vec;
- double *ros_vector;
- double **rcf_vector;
- complex<double> ***dc0m;
-
- int max_ici = 0;
+ int _idfc, _nxi;
+ int *_nshl_vec;
+ double *_xi_vec;
+ double *_ros_vec;
+ double **_rcf_vec;
+ complex<double> ***_dc0m;
+
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 configuration_count = 0;
- for (int i115 = 1; i115 <= nsph; i115++) {
- if (iog[i115 - 1] < i115) continue;
- configuration_count++;
- }
- 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[configuration_count]();
- ros_vector = new double[configuration_count]();
- rcf_vector = new double*[configuration_count];
- 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)); // was not from IOG
- input.read(reinterpret_cast<char *>(&(ros_vector[i115 - 1])), sizeof(double)); // was not from IOG
- int nsh = nshl_vector[i115 - 1]; // was not from IOG
- if (max_ici < (nsh + 1) / 2) max_ici = (nsh + 1) / 2;
- try {
- rcf_vector[i115 - 1] = new double[nsh](); // was not from IOG
- } 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)); // was not from IOG
- }
- }
- 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[iog[i162 - 1] - 1]; // was not from IOG
- 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.read(reinterpret_cast<char *>(&_nsph), sizeof(int));
+ input.read(reinterpret_cast<char *>(&_ies), sizeof(int));
+ _iog_vec = new int[_nsph]();
+ for (int i = 0; i < _nsph; i++)
+ input.read(reinterpret_cast<char *>(&(_iog_vec[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));
+ _xi_vec = new double[_nxi]();
+ for (int i = 0; i < _nxi; i++)
+ input.read(reinterpret_cast<char *>(&(_xi_vec[i])), sizeof(double));
+ int configurations = 0;
+ for (int i = 1; i <= _nsph; i++) {
+ if (_iog_vec[i - 1] >= i) configurations++;
+ }
+ _nshl_vec = new int[configurations]();
+ _ros_vec = new double[configurations]();
+ _rcf_vec = new double*[configurations];
+ for (int i115 = 1; i115 <= _nsph; i115++) {
+ if (_iog_vec[i115 - 1] < i115) continue;
+ input.read(reinterpret_cast<char *>(&(_nshl_vec[i115 - 1])), sizeof(int));
+ input.read(reinterpret_cast<char *>(&(_ros_vec[i115 - 1])), sizeof(double));
+ int nsh = _nshl_vec[i115 - 1];
+ if (i115 == 1) nsh += _ies;
+ _rcf_vec[i115 - 1] = new double[nsh]();
+ for (int nsi = 0; nsi < nsh; nsi++)
+ input.read(reinterpret_cast<char *>(&(_rcf_vec[i115 - 1][nsi])), sizeof(double));
+ }
+ _dc0m = new complex<double>**[configurations];
+ for (int di = 0; di < configurations; di++) {
+ _dc0m[di] = new complex<double>*[_nsph];
+ for (int dj = 0; dj < _nsph; dj++) _dc0m[di][dj] = new complex<double>[_nxi]();
+ } // di loop
+ for (int jxi468 = 1; jxi468 <= _nxi; jxi468++) {
+ if (_idfc != 0 && jxi468 > 1) continue;
+ for (int i162 = 1; i162 <= _nsph; 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 rval, ival;
+ input.read(reinterpret_cast<char *>(&rval), sizeof(double));
+ input.read(reinterpret_cast<char *>(&ival), sizeof(double));
+ _dc0m[i157][i162 - 1][jxi468 - 1] = complex<double>(rval, ival);
}
}
- input.close();
- } else { // Opening of the input file did not succeed.
- OpenConfigurationFileException ex(file_name);
- throw ex;
}
+ input.close();
ScattererConfiguration *conf = new ScattererConfiguration(
- nsph,
- xi_vec,
- nxi,
+ _nsph,
+ _xi_vec,
+ _nxi,
"XIV",
- iog,
- ros_vector,
- nshl_vector,
- rcf_vector,
+ _iog_vec,
+ _ros_vec,
+ _nshl_vec,
+ _rcf_vec,
_idfc,
- dc0m,
- false,
+ _dc0m,
+ (_ies == 1),
_exdc,
_wp,
_xip
@@ -889,12 +862,11 @@ void ScattererConfiguration::write_hdf5(string file_name) {
}
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));
+ output.write(reinterpret_cast<char *>(&ies), 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));