diff --git a/src/libnptm/Configuration.cpp b/src/libnptm/Configuration.cpp
index 517e7ea9726a80a758d512df4b69d1412bcb225f..7307de8c6e984f063b65c04dd334579fd508b0b8 100644
--- a/src/libnptm/Configuration.cpp
+++ b/src/libnptm/Configuration.cpp
@@ -13,921 +13,907 @@
using namespace std;
GeometryConfiguration::GeometryConfiguration(
- int _nsph, int _lm, int _in_pol, int _npnt, int _npntts, int _isam,
- int _li, int _le, int _mxndm, int _iavm,
- double *x, double *y, double *z,
- double in_th_start, double in_th_step, double in_th_end,
- double sc_th_start, double sc_th_step, double sc_th_end,
- double in_ph_start, double in_ph_step, double in_ph_end,
- double sc_ph_start, double sc_ph_step, double sc_ph_end,
- int _jwtm
- ) {
- number_of_spheres = _nsph;
- l_max = _lm;
- in_pol = _in_pol;
- npnt = _npnt;
- npntts = _npntts;
- meridional_type = _isam;
- li = _li;
- le = _le;
- mxndm = _mxndm;
- iavm = _iavm;
- in_theta_start = in_th_start;
- in_theta_step = in_th_step;
- in_theta_end = in_th_end;
- in_phi_start = in_ph_start;
- in_phi_step = in_ph_step;
- in_phi_end = in_ph_end;
- sc_theta_start = sc_th_start;
- sc_theta_step = sc_th_step;
- sc_theta_end = sc_th_end;
- sc_phi_start = sc_ph_start;
- sc_phi_step = sc_ph_step;
- sc_phi_end = sc_ph_end;
- jwtm = _jwtm;
- sph_x = x;
- sph_y = y;
- sph_z = z;
+ int _nsph, int _lm, int _in_pol, int _npnt, int _npntts, int _isam,
+ int _li, int _le, int _mxndm, int _iavm,
+ double *x, double *y, double *z,
+ double in_th_start, double in_th_step, double in_th_end,
+ double sc_th_start, double sc_th_step, double sc_th_end,
+ double in_ph_start, double in_ph_step, double in_ph_end,
+ double sc_ph_start, double sc_ph_step, double sc_ph_end,
+ int _jwtm
+ ) {
+ number_of_spheres = _nsph;
+ l_max = _lm;
+ in_pol = _in_pol;
+ npnt = _npnt;
+ npntts = _npntts;
+ meridional_type = _isam;
+ li = _li;
+ le = _le;
+ mxndm = _mxndm;
+ iavm = _iavm;
+ in_theta_start = in_th_start;
+ in_theta_step = in_th_step;
+ in_theta_end = in_th_end;
+ in_phi_start = in_ph_start;
+ in_phi_step = in_ph_step;
+ in_phi_end = in_ph_end;
+ sc_theta_start = sc_th_start;
+ sc_theta_step = sc_th_step;
+ sc_theta_end = sc_th_end;
+ sc_phi_start = sc_ph_start;
+ sc_phi_step = sc_ph_step;
+ sc_phi_end = sc_ph_end;
+ jwtm = _jwtm;
+ sph_x = x;
+ sph_y = y;
+ sph_z = z;
}
GeometryConfiguration::~GeometryConfiguration() {
- delete[] sph_x;
- delete[] sph_y;
- delete[] sph_z;
+ delete[] sph_x;
+ delete[] sph_y;
+ delete[] sph_z;
}
GeometryConfiguration* GeometryConfiguration::from_legacy(string file_name) {
- int num_lines = 0;
- int last_read_line = 0;
- string *file_lines;
- string str_target, str_num;
- smatch m;
- try {
- file_lines = load_file(file_name, &num_lines);
- } catch (exception &ex) {
- throw OpenConfigurationFileException(file_name);
- }
- int _nsph = 0, _lm = 0, _in_pol = 0, _npnt = 0, _npntts = 0, _isam = 0;
- int _li = 0, _le = 0, _mxndm = 0, _iavm = 0;
- regex re = regex("-?[0-9]+");
- str_target = file_lines[last_read_line++];
+ int num_lines = 0;
+ int last_read_line = 0;
+ string *file_lines;
+ string str_target, str_num;
+ smatch m;
+ try {
+ file_lines = load_file(file_name, &num_lines);
+ } catch (exception &ex) {
+ throw OpenConfigurationFileException(file_name);
+ }
+ int _nsph = 0, _lm = 0, _in_pol = 0, _npnt = 0, _npntts = 0, _isam = 0;
+ int _li = 0, _le = 0, _mxndm = 0, _iavm = 0;
+ regex re = regex("-?[0-9]+");
+ str_target = file_lines[last_read_line++];
+ regex_search(str_target, m, re);
+ _nsph = stoi(m.str());
+ if (_nsph == 1) {
+ for (int ri = 0; ri < 5; ri++) {
+ str_target = m.suffix().str();
+ regex_search(str_target, m, re);
+ if (ri == 0) _lm = stoi(m.str());
+ if (ri == 1) _in_pol = stoi(m.str());
+ if (ri == 2) _npnt = stoi(m.str());
+ if (ri == 3) _npntts = stoi(m.str());
+ if (ri == 4) _isam = stoi(m.str());
+ }
+ } else {
+ for (int ri = 0; ri < 8; ri++) {
+ str_target = m.suffix().str();
+ regex_search(str_target, m, re);
+ if (ri == 0) _li = stoi(m.str());
+ if (ri == 1) _le = stoi(m.str());
+ if (ri == 2) _mxndm = stoi(m.str());
+ if (ri == 3) _in_pol = stoi(m.str());
+ if (ri == 4) _npnt = stoi(m.str());
+ if (ri == 5) _npntts = stoi(m.str());
+ if (ri == 6) _iavm = stoi(m.str());
+ if (ri == 7) _isam = stoi(m.str());
+ }
+ }
+ double *x, *y, *z;
+ x = new double[_nsph];
+ y = new double[_nsph];
+ z = new double[_nsph];
+ if (_nsph == 1) {
+ x[0] = 0.0;
+ y[0] = 0.0;
+ z[0] = 0.0;
+ } else {
+ for (int i = 0; i < _nsph; i++) {
+ str_target = file_lines[last_read_line++];
+ re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
+ for (int ri = 0; ri < 3; ri++) {
regex_search(str_target, m, re);
- //sscanf(file_lines[last_read_line].c_str(), " %d", &_nsph);
- _nsph = stoi(m.str());
- if (_nsph == 1) {
- //sscanf(
- // file_lines[last_read_line++].c_str(),
- // " %*d %d %d %d %d %d",
- // &_lm, &_in_pol, &_npnt, &_npntts, &_isam
- //);
- for (int ri = 0; ri < 5; ri++) {
- str_target = m.suffix().str();
- regex_search(str_target, m, re);
- if (ri == 0) _lm = stoi(m.str());
- if (ri == 1) _in_pol = stoi(m.str());
- if (ri == 2) _npnt = stoi(m.str());
- if (ri == 3) _npntts = stoi(m.str());
- if (ri == 4) _isam = stoi(m.str());
- }
- } else {
- //sscanf(
- // file_lines[last_read_line++].c_str(),
- // " %*d %d %d %d %d %d %d %d %d",
- // &_li, &_le, &_mxndm, &_in_pol, &_npnt, &_npntts, &_iavm, &_isam
- //);
- for (int ri = 0; ri < 8; ri++) {
- str_target = m.suffix().str();
- regex_search(str_target, m, re);
- if (ri == 0) _li = stoi(m.str());
- if (ri == 1) _le = stoi(m.str());
- if (ri == 2) _mxndm = stoi(m.str());
- if (ri == 3) _in_pol = stoi(m.str());
- if (ri == 4) _npnt = stoi(m.str());
- if (ri == 5) _npntts = stoi(m.str());
- if (ri == 6) _iavm = stoi(m.str());
- if (ri == 7) _isam = stoi(m.str());
- }
- }
- double *x, *y, *z;
- x = new double[_nsph];
- y = new double[_nsph];
- z = new double[_nsph];
- if (_nsph == 1) {
- x[0] = 0.0;
- y[0] = 0.0;
- z[0] = 0.0;
- } else {
- for (int i = 0; i < _nsph; i++) {
- str_target = file_lines[last_read_line++];
- re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
- for (int ri = 0; ri < 3; ri++) {
- regex_search(str_target, m, re);
- str_num = regex_replace(m.str(), regex("D"), "e");
- str_num = regex_replace(str_num, regex("d"), "e");
- if (ri == 0) x[i] = stod(str_num);
- if (ri == 1) y[i] = stod(str_num);
- if (ri == 2) z[i] = stod(str_num);
- str_target = m.suffix().str();
- }
- }
- }
- double in_th_start, in_th_end, in_th_step, sc_th_start, sc_th_end, sc_th_step;
- re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
- str_target = file_lines[last_read_line++];
- for (int ri = 0; ri < 6; ri++) {
- regex_search(str_target, m, re);
- str_num = regex_replace(m.str(), regex("D"), "e");
- str_num = regex_replace(str_num, regex("d"), "e");
- if (ri == 0) in_th_start = stod(str_num);
- if (ri == 1) in_th_step = stod(str_num);
- if (ri == 2) in_th_end = stod(str_num);
- if (ri == 3) sc_th_start = stod(str_num);
- if (ri == 4) sc_th_step = stod(str_num);
- if (ri == 5) sc_th_end = stod(str_num);
- str_target = m.suffix().str();
- }
+ str_num = regex_replace(m.str(), regex("D"), "e");
+ str_num = regex_replace(str_num, regex("d"), "e");
+ if (ri == 0) x[i] = stod(str_num);
+ if (ri == 1) y[i] = stod(str_num);
+ if (ri == 2) z[i] = stod(str_num);
+ str_target = m.suffix().str();
+ }
+ }
+ }
+ double in_th_start, in_th_end, in_th_step, sc_th_start, sc_th_end, sc_th_step;
+ re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
+ str_target = file_lines[last_read_line++];
+ for (int ri = 0; ri < 6; ri++) {
+ regex_search(str_target, m, re);
+ str_num = regex_replace(m.str(), regex("D"), "e");
+ str_num = regex_replace(str_num, regex("d"), "e");
+ if (ri == 0) in_th_start = stod(str_num);
+ if (ri == 1) in_th_step = stod(str_num);
+ if (ri == 2) in_th_end = stod(str_num);
+ if (ri == 3) sc_th_start = stod(str_num);
+ if (ri == 4) sc_th_step = stod(str_num);
+ if (ri == 5) sc_th_end = stod(str_num);
+ str_target = m.suffix().str();
+ }
- double in_ph_start, in_ph_end, in_ph_step, sc_ph_start, sc_ph_end, sc_ph_step;
- str_target = file_lines[last_read_line++];
- for (int ri = 0; ri < 6; ri++) {
- regex_search(str_target, m, re);
- str_num = regex_replace(m.str(), regex("D"), "e");
- str_num = regex_replace(str_num, regex("d"), "e");
- if (ri == 0) in_ph_start = stod(str_num);
- if (ri == 1) in_ph_step = stod(str_num);
- if (ri == 2) in_ph_end = stod(str_num);
- if (ri == 3) sc_ph_start = stod(str_num);
- if (ri == 4) sc_ph_step = stod(str_num);
- if (ri == 5) sc_ph_end = stod(str_num);
- str_target = m.suffix().str();
- }
+ double in_ph_start, in_ph_end, in_ph_step, sc_ph_start, sc_ph_end, sc_ph_step;
+ str_target = file_lines[last_read_line++];
+ for (int ri = 0; ri < 6; ri++) {
+ regex_search(str_target, m, re);
+ str_num = regex_replace(m.str(), regex("D"), "e");
+ str_num = regex_replace(str_num, regex("d"), "e");
+ if (ri == 0) in_ph_start = stod(str_num);
+ if (ri == 1) in_ph_step = stod(str_num);
+ if (ri == 2) in_ph_end = stod(str_num);
+ if (ri == 3) sc_ph_start = stod(str_num);
+ if (ri == 4) sc_ph_step = stod(str_num);
+ if (ri == 5) sc_ph_end = stod(str_num);
+ str_target = m.suffix().str();
+ }
- int _jwtm;
- re = regex("[0-9]+");
- str_target = file_lines[last_read_line++];
- regex_search(str_target, m, re);
- //sscanf(file_lines[last_read_line].c_str(), " %d", &_nsph);
- _jwtm = stoi(m.str());
- GeometryConfiguration *conf = new GeometryConfiguration(
- _nsph, _lm, _in_pol, _npnt, _npntts, _isam,
- _li, _le, _mxndm, _iavm,
- x, y, z,
- in_th_start, in_th_step, in_th_end,
- sc_th_start, sc_th_step, sc_th_end,
- in_ph_start, in_ph_step, in_ph_end,
- sc_ph_start, sc_ph_step, sc_ph_end,
- _jwtm
- );
- return conf;
+ int _jwtm;
+ re = regex("[0-9]+");
+ str_target = file_lines[last_read_line++];
+ regex_search(str_target, m, re);
+ _jwtm = stoi(m.str());
+ GeometryConfiguration *conf = new GeometryConfiguration(
+ _nsph, _lm, _in_pol, _npnt, _npntts, _isam,
+ _li, _le, _mxndm, _iavm,
+ x, y, z,
+ in_th_start, in_th_step, in_th_end,
+ sc_th_start, sc_th_step, sc_th_end,
+ in_ph_start, in_ph_step, in_ph_end,
+ sc_ph_start, sc_ph_step, sc_ph_end,
+ _jwtm
+ );
+ return conf;
}
ScattererConfiguration::ScattererConfiguration(
- int nsph,
- double *scale_vector,
- int nxi,
- string variable_name,
- int *iog_vector,
- double *ros_vector,
- int *nshl_vector,
- double **rcf_vector,
- int dielectric_func_type,
- complex<double> ***dc_matrix,
- bool is_external,
- double ex,
- double w,
- double x
- ) {
- number_of_spheres = nsph;
- number_of_scales = nxi;
- reference_variable_name = variable_name;
- iog_vec = iog_vector;
- radii_of_spheres = ros_vector;
- nshl_vec = nshl_vector;
- rcf = rcf_vector;
- idfc = dielectric_func_type,
- dc0_matrix = dc_matrix;
- use_external_sphere = is_external;
- exdc = ex;
- wp = w;
- xip = x;
- if (variable_name == "XIV") scale_vec = scale_vector;
- else {
- scale_vec = new double[number_of_scales]();
- const double pi2 = 2.0 * acos(-1.0);
- const double evc = 6.5821188e-16;
- for (int si = 0; si < number_of_scales; si++) {
- double value = scale_vector[si];
- if (variable_name.compare("WNS") == 0) value *= (3.0e8 / wp);
- else if (variable_name.compare("WLS") == 0) value = pi2 / value * 3.0e8 / wp;
- else if (variable_name.compare("PUS") == 0) value /= wp;
- else if (variable_name.compare("EVS") == 0) value /= (evc * wp);
- scale_vec[si] = value;
- }
- }
+ int nsph,
+ double *scale_vector,
+ int nxi,
+ string variable_name,
+ int *iog_vector,
+ double *ros_vector,
+ int *nshl_vector,
+ double **rcf_vector,
+ int dielectric_func_type,
+ complex<double> ***dc_matrix,
+ bool is_external,
+ double ex,
+ double w,
+ double x
+ ) {
+ number_of_spheres = nsph;
+ number_of_scales = nxi;
+ reference_variable_name = variable_name;
+ iog_vec = iog_vector;
+ radii_of_spheres = ros_vector;
+ nshl_vec = nshl_vector;
+ rcf = rcf_vector;
+ idfc = dielectric_func_type,
+ dc0_matrix = dc_matrix;
+ use_external_sphere = is_external;
+ exdc = ex;
+ wp = w;
+ xip = x;
+ if (variable_name == "XIV") scale_vec = scale_vector;
+ else {
+ scale_vec = new double[number_of_scales]();
+ const double pi2 = 2.0 * acos(-1.0);
+ const double evc = 6.5821188e-16;
+ for (int si = 0; si < number_of_scales; si++) {
+ double value = scale_vector[si];
+ if (variable_name.compare("WNS") == 0) value *= (3.0e8 / wp);
+ else if (variable_name.compare("WLS") == 0) value = pi2 / value * 3.0e8 / wp;
+ else if (variable_name.compare("PUS") == 0) value /= wp;
+ else if (variable_name.compare("EVS") == 0) value /= (evc * wp);
+ scale_vec[si] = value;
+ }
+ }
}
ScattererConfiguration::~ScattererConfiguration() {
- int max_ici = 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;
- }
- for (int i = 0; i < max_ici; 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++) {
- delete[] rcf[i];
- }
- delete[] nshl_vec;
- delete[] radii_of_spheres;
- delete[] iog_vec;
- delete[] scale_vec;
+ int max_ici = 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;
+ }
+ for (int i = 0; i < max_ici; 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++) {
+ delete[] rcf[i];
+ }
+ delete[] nshl_vec;
+ delete[] radii_of_spheres;
+ delete[] iog_vec;
+ delete[] scale_vec;
}
ScattererConfiguration* ScattererConfiguration::from_binary(string file_name, string mode) {
- int nsph;
- int *iog;
- double _exdc, _wp, _xip;
- int _idfc, nxi;
- int *nshl_vector;
- double *xi_vec;
- 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));
- }
- 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 (bad_alloc &ex) {
- throw UnrecognizedConfigurationException("Wrong parameter set: invalid number of scales " + 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 *>(&(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 (bad_alloc &ex) {
- throw UnrecognizedConfigurationException("Wrong parameter set: invalid number of layers " + 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]();
- }
- }
- 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.
- throw OpenConfigurationFileException(file_name);
- }
- } else { // A different binary format was chosen.
- //TODO: this part is not yet implemented.
- // Functions to write optimized file formats may be invoked here.
+ int nsph;
+ int *iog;
+ double _exdc, _wp, _xip;
+ int _idfc, nxi;
+ int *nshl_vector;
+ double *xi_vec;
+ 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));
+ }
+ 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 (bad_alloc &ex) {
+ throw UnrecognizedConfigurationException("Wrong parameter set: invalid number of scales " + 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;
}
- ScattererConfiguration *conf = new ScattererConfiguration(
- nsph,
- xi_vec,
- nxi,
- "XIV",
- iog,
- ros_vector,
- nshl_vector,
- rcf_vector,
- _idfc,
- dc0m,
- false,
- _exdc,
- _wp,
- _xip
- );
- return conf;
-}
-
-ScattererConfiguration* ScattererConfiguration::from_dedfb(string dedfb_file_name) {
- int num_lines = 0;
- int last_read_line = 0;
- string *file_lines;
- regex re;
- smatch m;
+ 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 {
- file_lines = load_file(dedfb_file_name, &num_lines);
- } catch (exception &ex) {
- throw OpenConfigurationFileException(dedfb_file_name);
+ rcf_vector[i115 - 1] = new double[nsh]();
+ } catch (bad_alloc &ex) {
+ throw UnrecognizedConfigurationException("Wrong parameter set: invalid number of layers " + nsh);
}
- int nsph, ies;
- int max_ici = 0;
- re = regex("[0-9]+");
- string str_target = file_lines[last_read_line];
- for (int ri = 0; ri < 2; ri++) {
- regex_search(str_target, m, re);
- if (ri == 0) nsph = stoi(m.str());
- if (ri == 1) ies = stoi(m.str());
- str_target = m.suffix().str();
+ for (int nsi = 0; nsi < nsh; nsi++) {
+ input.read(reinterpret_cast<char *>(&(rcf_vector[i115 - 1][nsi])), sizeof(double));
}
- if (ies != 0) ies = 1;
- double _exdc, _wp, _xip;
- str_target = file_lines[++last_read_line];
- re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
- for (int ri = 0; ri < 3; ri++) {
- 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");
- if (ri == 0) _exdc = stod(str_number);
- if (ri == 1) _wp = stod(str_number);
- if (ri == 2) _xip = stod(str_number);
- str_target = m.suffix().str();
+ }
+ 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]();
}
- int _idfc, nxi, instpc, insn;
- re = regex("-?[0-9]+");
- for (int ri = 0; ri < 4; ri++) {
- regex_search(str_target, m, re);
- if (ri == 0) _idfc = stoi(m.str());
- if (ri == 1) nxi = stoi(m.str());
- if (ri == 2) instpc = stoi(m.str());
- if (ri == 3) insn = stoi(m.str());
- str_target = m.suffix().str();
+ }
+ 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.
+ throw OpenConfigurationFileException(file_name);
+ }
+ } else { // A different binary format was chosen.
+ //TODO: this part is not yet implemented.
+ // Functions to write optimized file formats may be invoked here.
+ }
+ ScattererConfiguration *conf = new ScattererConfiguration(
+ nsph,
+ xi_vec,
+ nxi,
+ "XIV",
+ iog,
+ ros_vector,
+ nshl_vector,
+ rcf_vector,
+ _idfc,
+ dc0m,
+ false,
+ _exdc,
+ _wp,
+ _xip
+ );
+ return conf;
+}
+
+ScattererConfiguration* ScattererConfiguration::from_dedfb(string dedfb_file_name) {
+ int num_lines = 0;
+ int last_read_line = 0;
+ string *file_lines;
+ regex re;
+ smatch m;
+ try {
+ file_lines = load_file(dedfb_file_name, &num_lines);
+ } catch (exception &ex) {
+ throw OpenConfigurationFileException(dedfb_file_name);
+ }
+ int nsph, ies;
+ int max_ici = 0;
+ re = regex("[0-9]+");
+ string str_target = file_lines[last_read_line];
+ for (int ri = 0; ri < 2; ri++) {
+ regex_search(str_target, m, re);
+ if (ri == 0) nsph = stoi(m.str());
+ if (ri == 1) ies = stoi(m.str());
+ str_target = m.suffix().str();
+ }
+ if (ies != 0) ies = 1;
+ double _exdc, _wp, _xip;
+ str_target = file_lines[++last_read_line];
+ re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
+ for (int ri = 0; ri < 3; ri++) {
+ 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");
+ if (ri == 0) _exdc = stod(str_number);
+ if (ri == 1) _wp = stod(str_number);
+ if (ri == 2) _xip = stod(str_number);
+ str_target = m.suffix().str();
+ }
+ int _idfc, nxi, instpc, insn;
+ re = regex("-?[0-9]+");
+ for (int ri = 0; ri < 4; ri++) {
+ regex_search(str_target, m, re);
+ if (ri == 0) _idfc = stoi(m.str());
+ if (ri == 1) nxi = stoi(m.str());
+ if (ri == 2) instpc = stoi(m.str());
+ if (ri == 3) insn = stoi(m.str());
+ str_target = m.suffix().str();
+ }
- double *variable_vector;
- string variable_name;
+ double *variable_vector;
+ string variable_name;
- if (_idfc < 0) { // Diel. functions at XIP value and XI is scale factor
- variable_name = "XIV";
- if (instpc < 1) { // The variable vector is explicitly defined.
- double xi;
- //int xi_exp;
- List<double> xi_vector;
- //sscanf(file_lines[++last_read_line].c_str(), " %9lE D%d", &xi, &xi_exp);
- str_target = file_lines[++last_read_line];
- 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");
- xi = stod(str_number);
- xi_vector.set(0, xi);
- for (int jxi310 = 1; jxi310 < nxi; jxi310++) {
- 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");
- xi = stod(str_number);
- xi_vector.append(xi);
- }
- variable_vector = xi_vector.to_array();
- } else { // instpc >= 1: the variable vector is defined in steps
- double xi, xi_step;
- str_target = file_lines[++last_read_line];
- re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
- regex_search(str_target, m, re);
- for (int ri = 0; ri < 2; ri++) {
- 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");
- if (ri == 0) xi = stod(str_number);
- if (ri == 1) xi_step = stod(str_number);
- str_target = m.suffix().str();
- }
- variable_vector = new double[nxi]();
- for (int jxi320 = 0; jxi320 < nxi; jxi320++) {
- variable_vector[jxi320] = xi;
- xi += xi_step;
- }
- }
- } else { // idfc >= 0
- variable_vector = new double[nxi]();
- if (instpc == 0) { // The variable vector is explicitly defined
- double vs;
- //int vs_exp;
- for (int jxi_r = 0; jxi_r < nxi; jxi_r++) {
- //sscanf(file_lines[++last_read_line].c_str(), " %lf D%d", &vs, &vs_exp);
- str_target = file_lines[++last_read_line];
- 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");
- vs = stod(str_number);
- variable_vector[jxi_r] = vs;
- }
- switch (insn) {
- case 1: //xi vector definition
- variable_name = "XIV";
- break;
- case 2: //wave number vector definition
- variable_name = "WNS";
- break;
- case 3: //wavelength vector definition
- variable_name = "WLS";
- break;
- case 4: //pu vector definition
- variable_name = "PUS";
- break;
- case 5: //eV vector definition
- variable_name = "EVS";
- break;
- }
- } else { // The variable vector needs to be computed in steps
- double vs, vs_step;
- //int vs_exp, vs_step_exp;
- //sscanf(file_lines[++last_read_line].c_str(), " %lf D%d %lf D%d", &vs, &vs_exp, &vs_step, &vs_step_exp);
- str_target = file_lines[++last_read_line];
- re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
- regex_search(str_target, m, re);
- for (int ri = 0; ri < 2; ri++) {
- 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");
- if (ri == 0) vs = stod(str_number);
- if (ri == 1) vs_step = stod(str_number);
- str_target = m.suffix().str();
- }
- for (int jxi110w = 0; jxi110w < nxi; jxi110w++) {
- variable_vector[jxi110w] = vs;
- vs += vs_step;
- }
- switch (insn) {
- case 1: //xi vector definition
- variable_name = "XIV";
- break;
- case 2: //wave number vector definition
- variable_name = "WNS";
- break;
- case 3: //wavelength vector definition
- variable_name = "WLS";
- break;
- case 4: //pu vector definition
- variable_name = "PUS";
- break;
- case 5: //eV vector definition
- variable_name = "EVS";
- break;
- }
- }
- }
- 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]();
- //printf("\nDEBUG: reading IOG from %s", file_lines[last_read_line].c_str());
- //fflush(stdout);
- 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++;
- //printf("DEBUG: reading IOG from %s", file_lines[last_read_line].c_str());
- //fflush(stdout);
- }
- }
- for (int i113 = 1; i113 <= nsph; i113++) {
- int i_val, nsh;
- double ros_val;
- //int ros_val_exp;
- if (iog_vector[i113 - 1] < i113) {
- rcf_vector[i113 - 1] = new double[1]();
- continue;
- }
- //sscanf(file_lines[++last_read_line].c_str(), " %d %lf D%d", &i_val, &ros_val, &ros_val_exp);
- re = regex("[0-9]+");
- str_target = file_lines[++last_read_line];
- regex_search(str_target, m, re);
- i_val = stoi(m.str());
- str_target = m.suffix();
- 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];
- if (i113 == 1) nsh += ies;
- 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();
- 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];
- 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]();
- }
- }
+ if (_idfc < 0) { // Diel. functions at XIP value and XI is scale factor
+ variable_name = "XIV";
+ if (instpc < 1) { // The variable vector is explicitly defined.
+ double xi;
+ List<double> xi_vector;
+ str_target = file_lines[++last_read_line];
+ 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");
+ xi = stod(str_number);
+ xi_vector.set(0, xi);
+ for (int jxi310 = 1; jxi310 < nxi; jxi310++) {
+ 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");
+ xi = stod(str_number);
+ xi_vector.append(xi);
+ }
+ variable_vector = xi_vector.to_array();
+ } else { // instpc >= 1: the variable vector is defined in steps
+ double xi, xi_step;
+ str_target = file_lines[++last_read_line];
+ re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
+ regex_search(str_target, m, re);
+ for (int ri = 0; ri < 2; ri++) {
+ 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");
+ if (ri == 0) xi = stod(str_number);
+ if (ri == 1) xi_step = stod(str_number);
+ str_target = m.suffix().str();
+ }
+ variable_vector = new double[nxi]();
+ for (int jxi320 = 0; jxi320 < nxi; jxi320++) {
+ variable_vector[jxi320] = xi;
+ xi += xi_step;
+ }
+ }
+ } else { // idfc >= 0
+ variable_vector = new double[nxi]();
+ if (instpc == 0) { // The variable vector is explicitly defined
+ double vs;
+ for (int jxi_r = 0; jxi_r < nxi; jxi_r++) {
+ str_target = file_lines[++last_read_line];
re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
- for (int jxi468 = 1; jxi468 <= nxi; jxi468++) {
- if (_idfc != 0 && jxi468 > 1) continue;
- 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;
- //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();
- 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] = std::complex<double>(dc0_real, dc0_img);
- }
- }
- }
+ 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");
+ vs = stod(str_number);
+ variable_vector[jxi_r] = vs;
+ }
+ switch (insn) {
+ case 1: //xi vector definition
+ variable_name = "XIV";
+ break;
+ case 2: //wave number vector definition
+ variable_name = "WNS";
+ break;
+ case 3: //wavelength vector definition
+ variable_name = "WLS";
+ break;
+ case 4: //pu vector definition
+ variable_name = "PUS";
+ break;
+ case 5: //eV vector definition
+ variable_name = "EVS";
+ break;
+ }
+ } else { // The variable vector needs to be computed in steps
+ double vs, vs_step;
+ str_target = file_lines[++last_read_line];
+ re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
+ regex_search(str_target, m, re);
+ for (int ri = 0; ri < 2; ri++) {
+ 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");
+ if (ri == 0) vs = stod(str_number);
+ if (ri == 1) vs_step = stod(str_number);
+ str_target = m.suffix().str();
+ }
+ for (int jxi110w = 0; jxi110w < nxi; jxi110w++) {
+ variable_vector[jxi110w] = vs;
+ vs += vs_step;
+ }
+ switch (insn) {
+ case 1: //xi vector definition
+ variable_name = "XIV";
+ break;
+ case 2: //wave number vector definition
+ variable_name = "WNS";
+ break;
+ case 3: //wavelength vector definition
+ variable_name = "WLS";
+ break;
+ case 4: //pu vector definition
+ variable_name = "PUS";
+ break;
+ case 5: //eV vector definition
+ variable_name = "EVS";
+ break;
+ }
+ }
+ }
+ //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) {
+ string str_target = file_lines[++last_read_line];
+ while(regex_search(str_target, m, re)) {
+ iog_vector[filled_iogs++] = stoi(m.str());
+ str_target = m.suffix().str();
+ if (filled_iogs == nsph) break;
+ }
+ }
+ 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]+");
+ str_target = file_lines[++last_read_line];
+ regex_search(str_target, m, re);
+ i_val = stoi(m.str());
+ str_target = m.suffix();
+ 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];
+ if (i113 == 1) nsh += ies;
+ 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();
+ 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];
+ 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]();
+ }
+ }
+ re = regex("-?[0-9]+\\.[0-9]+([eEdD][-+]?)?[0-9]+");
+ for (int jxi468 = 1; jxi468 <= nxi; jxi468++) {
+ if (_idfc != 0 && jxi468 > 1) continue;
+ 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;
+ //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();
+ 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] = std::complex<double>(dc0_real, dc0_img);
+ }
+ }
+ }
- ScattererConfiguration *config = new ScattererConfiguration(
- nsph,
- variable_vector,
- nxi,
- variable_name,
- iog_vector,
- ros_vector,
- nshl_vector,
- rcf_vector,
- _idfc,
- dc0m,
- (ies > 0),
- _exdc,
- _wp,
- _xip
- );
- delete[] file_lines;
- return config;
+ ScattererConfiguration *config = new ScattererConfiguration(
+ nsph,
+ variable_vector,
+ nxi,
+ variable_name,
+ iog_vector,
+ ros_vector,
+ nshl_vector,
+ rcf_vector,
+ _idfc,
+ dc0m,
+ (ies > 0),
+ _exdc,
+ _wp,
+ _xip
+ );
+ delete[] file_lines;
+ return config;
}
void ScattererConfiguration::print() {
- int ies = (use_external_sphere)? 1 : 0;
- int max_ici = 0;
- printf("### CONFIGURATION DATA ###\n");
- printf("NSPH = %d\n", number_of_spheres);
- printf("ROS = [");
- for (int i = 0; i < number_of_spheres; i++) printf("\t%lg", radii_of_spheres[i]);
- printf("\t]\n");
- printf("IOG = [");
- for (int i = 0; i < number_of_spheres; i++) printf("\t%d", iog_vec[i]);
- printf("\t]\n");
- printf("NSHL = [");
- for (int i = 0; i < number_of_spheres; i++) printf("\t%d", nshl_vec[i]);
- printf("\t]\n");
- printf("RCF = [\n");
- for (int i = 1; i <= number_of_spheres; i++) {
- int nsh = nshl_vec[i - 1];
- if (i == 1) nsh += ies;
- if (max_ici < (nsh + 1) / 2) max_ici = (nsh + 1) / 2;
- printf(" [");
- for (int ns = 0; ns < nsh; ns++) {
- printf("\t%lg", rcf[i - 1][ns]);
- }
- printf("\t]\n");
- }
- printf(" ]\n");
- printf("SCALE = %s\n", reference_variable_name.c_str());
- printf("NXI = %d\n", number_of_scales);
- printf("VEC = [");
- for (int i = 0; i < number_of_scales; i++) printf("\t%lg", scale_vec[i]);
- printf("\t]\n");
- printf("DC0M = [\n");
- for (int i = 0; i < max_ici; i++) {
- printf(" [\n");
- for (int j = 0; j < number_of_spheres; j++) {
- printf(" [");
- for (int k = 0; k < number_of_scales; k++) {
- if (idfc != 0 and k > 0) continue;
- printf("\t%lg + i(%lg)", dc0_matrix[i][j][k].real(), dc0_matrix[i][j][k].imag());
- }
- printf("\t]\n");
- }
- printf(" ]\n");
- }
- printf(" ]\n");
+ int ies = (use_external_sphere)? 1 : 0;
+ int max_ici = 0;
+ printf("### CONFIGURATION DATA ###\n");
+ printf("NSPH = %d\n", number_of_spheres);
+ printf("ROS = [");
+ for (int i = 0; i < number_of_spheres; i++) printf("\t%lg", radii_of_spheres[i]);
+ printf("\t]\n");
+ printf("IOG = [");
+ for (int i = 0; i < number_of_spheres; i++) printf("\t%d", iog_vec[i]);
+ printf("\t]\n");
+ printf("NSHL = [");
+ for (int i = 0; i < number_of_spheres; i++) printf("\t%d", nshl_vec[i]);
+ printf("\t]\n");
+ printf("RCF = [\n");
+ for (int i = 1; i <= number_of_spheres; i++) {
+ int nsh = nshl_vec[i - 1];
+ if (i == 1) nsh += ies;
+ if (max_ici < (nsh + 1) / 2) max_ici = (nsh + 1) / 2;
+ printf(" [");
+ for (int ns = 0; ns < nsh; ns++) {
+ printf("\t%lg", rcf[i - 1][ns]);
+ }
+ printf("\t]\n");
+ }
+ printf(" ]\n");
+ printf("SCALE = %s\n", reference_variable_name.c_str());
+ printf("NXI = %d\n", number_of_scales);
+ printf("VEC = [");
+ for (int i = 0; i < number_of_scales; i++) printf("\t%lg", scale_vec[i]);
+ printf("\t]\n");
+ printf("DC0M = [\n");
+ for (int i = 0; i < max_ici; i++) {
+ printf(" [\n");
+ for (int j = 0; j < number_of_spheres; j++) {
+ printf(" [");
+ for (int k = 0; k < number_of_scales; k++) {
+ if (idfc != 0 and k > 0) continue;
+ printf("\t%lg + i(%lg)", dc0_matrix[i][j][k].real(), dc0_matrix[i][j][k].imag());
+ }
+ printf("\t]\n");
+ }
+ printf(" ]\n");
+ }
+ printf(" ]\n");
}
void ScattererConfiguration::write_binary(string file_name, string mode) {
- const double two_pi = acos(0.0) * 4.0;
- const double evc = 6.5821188e-16;
- 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;
- double *xi_vec;
- if (reference_variable_name.compare("XIV") == 0) xi_vec = scale_vec;
- else {
- is_new_vector = true;
- xi_vec = new double[number_of_scales];
- if (reference_variable_name.compare("WNS") == 0) {
- for (int i = 0; i < number_of_scales; i++)
- xi_vec[i] = 3.0e8 * scale_vec[i] / wp;
- } else if (reference_variable_name.compare("WLS") == 0) {
- for (int i = 0; i < number_of_scales; i++) {
- double wn = two_pi / scale_vec[i];
- xi_vec[i] = 3.0e8 * wn / wp;
- }
- } else if (reference_variable_name.compare("PUS") == 0) {
- for (int i = 0; i < number_of_scales; i++)
- xi_vec[i] = scale_vec[i] / wp;
- } else if (reference_variable_name.compare("EVS") == 0) {
- for (int i = 0; i < number_of_scales; i++) {
- double pu = scale_vec[i] / evc;
- xi_vec[i] = pu / wp;
- }
- } else {
- throw UnrecognizedConfigurationException(
- "Wrong parameter set: unrecognized scale type "
- + reference_variable_name
- );
- }
- }
- 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 *>(&(xi_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));
- }
- 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));
- }
- }
- }
- if (is_new_vector) delete[] xi_vec;
- output.close();
+ const double two_pi = acos(0.0) * 4.0;
+ const double evc = 6.5821188e-16;
+ 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;
+ double *xi_vec;
+ if (reference_variable_name.compare("XIV") == 0) xi_vec = scale_vec;
+ else {
+ is_new_vector = true;
+ xi_vec = new double[number_of_scales];
+ if (reference_variable_name.compare("WNS") == 0) {
+ for (int i = 0; i < number_of_scales; i++)
+ xi_vec[i] = 3.0e8 * scale_vec[i] / wp;
+ } else if (reference_variable_name.compare("WLS") == 0) {
+ for (int i = 0; i < number_of_scales; i++) {
+ double wn = two_pi / scale_vec[i];
+ xi_vec[i] = 3.0e8 * wn / wp;
+ }
+ } else if (reference_variable_name.compare("PUS") == 0) {
+ for (int i = 0; i < number_of_scales; i++)
+ xi_vec[i] = scale_vec[i] / wp;
+ } else if (reference_variable_name.compare("EVS") == 0) {
+ for (int i = 0; i < number_of_scales; i++) {
+ double pu = scale_vec[i] / evc;
+ xi_vec[i] = pu / wp;
}
+ } else {
+ throw UnrecognizedConfigurationException(
+ "Wrong parameter set: unrecognized scale type "
+ + reference_variable_name
+ );
+ }
+ }
+ 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 *>(&(xi_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));
+ }
+ 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));
+ }
+ }
+ }
+ if (is_new_vector) delete[] xi_vec;
+ output.close();
+ }
}
void ScattererConfiguration::write_formatted(string file_name) {
- const double evc = 6.5821188e-16;
- const double two_pi = acos(0.0) * 4.0;
- double *xi_vec;
- int ici;
- int ies = (use_external_sphere)? 1: 0;
- FILE *output = fopen(file_name.c_str(), "w");
- int scale_type = -1;
- if (reference_variable_name.compare("XIV") == 0) scale_type = 0;
- else if (reference_variable_name.compare("WNS") == 0) scale_type = 1;
- else if (reference_variable_name.compare("WLS") == 0) scale_type = 2;
- else if (reference_variable_name.compare("PUS") == 0) scale_type = 3;
- else if (reference_variable_name.compare("EVS") == 0) scale_type = 4;
- if (idfc >= 0) { // Dielectric functions are constant or they depend on XI
- double *pu_vec, *ev_vec, *wn_vec, *wl_vec;
- xi_vec = new double[number_of_scales];
- pu_vec = new double[number_of_scales];
- ev_vec = new double[number_of_scales];
- wn_vec = new double[number_of_scales];
- wl_vec = new double[number_of_scales];
- switch (scale_type) {
- case 0:
- fprintf(output, " JXI XIV WNS WLS PUS EVS\n");
- for (int i = 0; i < number_of_scales; i++) {
- xi_vec[i] = scale_vec[i];
- pu_vec[i] = xi_vec[i] * wp;
- ev_vec[i] = pu_vec[i] * evc;
- wn_vec[i] = pu_vec[i] / 3.0e8;
- wl_vec[i] = two_pi / wn_vec[i];
- fprintf(
- output,
- "%5d%13.4lE%13.4lE%13.4lE%13.4lE%13.4lE\n",
- (i + 1),
- xi_vec[i],
- wn_vec[i],
- wl_vec[i],
- pu_vec[i],
- ev_vec[i]
- );
- }
- break;
- case 1:
- fprintf(output, " JXI WNS WLS PUS EVS XIV\n");
- for (int i = 0; i < number_of_scales; i++) {
- wn_vec[i] = scale_vec[i];
- wl_vec[i] = two_pi / wn_vec[i];
- xi_vec[i] = 3.0e8 * wn_vec[i] / wp;
- pu_vec[i] = xi_vec[i] * wp;
- ev_vec[i] = pu_vec[i] * evc;
- fprintf(
- output,
- "%5d%13.4lE%13.4lE%13.4lE%13.4lE%13.4lE\n",
- (i + 1),
- wn_vec[i],
- wl_vec[i],
- pu_vec[i],
- ev_vec[i],
- xi_vec[i]
- );
- }
- break;
- case 2:
- fprintf(output, " JXI WLS WNS PUS EVS XIV\n");
- for (int i = 0; i < number_of_scales; i++) {
- wl_vec[i] = scale_vec[i];
- wn_vec[i] = two_pi / wl_vec[i];
- xi_vec[i] = 3.0e8 * wn_vec[i] / wp;
- pu_vec[i] = xi_vec[i] * wp;
- ev_vec[i] = pu_vec[i] * evc;
- fprintf(
- output,
- "%5d%13.4lE%13.4lE%13.4lE%13.4lE%13.4lE\n",
- (i + 1),
- wl_vec[i],
- wn_vec[i],
- pu_vec[i],
- ev_vec[i],
- xi_vec[i]
- );
- }
- break;
- case 3:
- fprintf(output, " JXI PUS WNS WLS EVS XIV\n");
- for (int i = 0; i < number_of_scales; i++) {
- pu_vec[i] = scale_vec[i];
- xi_vec[i] = pu_vec[i] / wp;
- wn_vec[i] = pu_vec[i] / 3.0e8;
- wl_vec[i] = two_pi / wn_vec[i];
- ev_vec[i] = pu_vec[i] * evc;
- fprintf(
- output,
- "%5d%13.4lE%13.4lE%13.4lE%13.4lE%13.4lE\n",
- (i + 1),
- pu_vec[i],
- wn_vec[i],
- wl_vec[i],
- ev_vec[i],
- xi_vec[i]
- );
- }
- break;
- case 4:
- fprintf(output, " JXI EVS WNS WLS PUS XIV\n");
- for (int i = 0; i < number_of_scales; i++) {
- ev_vec[i] = scale_vec[i];
- pu_vec[i] = ev_vec[i] / evc;
- xi_vec[i] = pu_vec[i] / wp;
- wn_vec[i] = pu_vec[i] / 3.0e8;
- wl_vec[i] = two_pi / wn_vec[i];
- fprintf(
- output,
- "%5d%13.4lE%13.4lE%13.4lE%13.4lE%13.4lE\n",
- (i + 1),
- ev_vec[i],
- wn_vec[i],
- wl_vec[i],
- pu_vec[i],
- xi_vec[i]
- );
- }
- break;
- default:
- throw UnrecognizedConfigurationException(
- "Wonrg parameter set: unrecognized scale variable type "
- + reference_variable_name
- );
- break;
- }
- } else { // idfc < 0, Dielectric functions are at XIP and XI is scale for dimensions
- double pu, wn;
- xi_vec = scale_vec;
- pu = xip * wp;
- wn = pu / 3.0e8;
- fprintf(output, " XIP WN WL PU EV\n");
- fprintf(output, " %13.4lE", xip);
- fprintf(output, "%13.4lE", wn);
- fprintf(output, "%13.4lE", two_pi / wn);
- fprintf(output, "%13.4lE", pu);
- fprintf(output, "%13.4lE\n", pu * evc);
- fprintf(output, " SCALE FACTORS XI\n");
- for (int i = 0; i < number_of_scales; i++)
- fprintf(output, "%5d%13.4lE\n", (i + 1), xi_vec[i]);
+ const double evc = 6.5821188e-16;
+ const double two_pi = acos(0.0) * 4.0;
+ double *xi_vec;
+ int ici;
+ int ies = (use_external_sphere)? 1: 0;
+ FILE *output = fopen(file_name.c_str(), "w");
+ int scale_type = -1;
+ if (reference_variable_name.compare("XIV") == 0) scale_type = 0;
+ else if (reference_variable_name.compare("WNS") == 0) scale_type = 1;
+ else if (reference_variable_name.compare("WLS") == 0) scale_type = 2;
+ else if (reference_variable_name.compare("PUS") == 0) scale_type = 3;
+ else if (reference_variable_name.compare("EVS") == 0) scale_type = 4;
+ if (idfc >= 0) { // Dielectric functions are constant or they depend on XI
+ double *pu_vec, *ev_vec, *wn_vec, *wl_vec;
+ xi_vec = new double[number_of_scales];
+ pu_vec = new double[number_of_scales];
+ ev_vec = new double[number_of_scales];
+ wn_vec = new double[number_of_scales];
+ wl_vec = new double[number_of_scales];
+ switch (scale_type) {
+ case 0:
+ fprintf(output, " JXI XIV WNS WLS PUS EVS\n");
+ for (int i = 0; i < number_of_scales; i++) {
+ xi_vec[i] = scale_vec[i];
+ pu_vec[i] = xi_vec[i] * wp;
+ ev_vec[i] = pu_vec[i] * evc;
+ wn_vec[i] = pu_vec[i] / 3.0e8;
+ wl_vec[i] = two_pi / wn_vec[i];
+ fprintf(
+ output,
+ "%5d%13.4lE%13.4lE%13.4lE%13.4lE%13.4lE\n",
+ (i + 1),
+ xi_vec[i],
+ wn_vec[i],
+ wl_vec[i],
+ pu_vec[i],
+ ev_vec[i]
+ );
+ }
+ break;
+ case 1:
+ fprintf(output, " JXI WNS WLS PUS EVS XIV\n");
+ for (int i = 0; i < number_of_scales; i++) {
+ wn_vec[i] = scale_vec[i];
+ wl_vec[i] = two_pi / wn_vec[i];
+ xi_vec[i] = 3.0e8 * wn_vec[i] / wp;
+ pu_vec[i] = xi_vec[i] * wp;
+ ev_vec[i] = pu_vec[i] * evc;
+ fprintf(
+ output,
+ "%5d%13.4lE%13.4lE%13.4lE%13.4lE%13.4lE\n",
+ (i + 1),
+ wn_vec[i],
+ wl_vec[i],
+ pu_vec[i],
+ ev_vec[i],
+ xi_vec[i]
+ );
+ }
+ break;
+ case 2:
+ fprintf(output, " JXI WLS WNS PUS EVS XIV\n");
+ for (int i = 0; i < number_of_scales; i++) {
+ wl_vec[i] = scale_vec[i];
+ wn_vec[i] = two_pi / wl_vec[i];
+ xi_vec[i] = 3.0e8 * wn_vec[i] / wp;
+ pu_vec[i] = xi_vec[i] * wp;
+ ev_vec[i] = pu_vec[i] * evc;
+ fprintf(
+ output,
+ "%5d%13.4lE%13.4lE%13.4lE%13.4lE%13.4lE\n",
+ (i + 1),
+ wl_vec[i],
+ wn_vec[i],
+ pu_vec[i],
+ ev_vec[i],
+ xi_vec[i]
+ );
+ }
+ break;
+ case 3:
+ fprintf(output, " JXI PUS WNS WLS EVS XIV\n");
+ for (int i = 0; i < number_of_scales; i++) {
+ pu_vec[i] = scale_vec[i];
+ xi_vec[i] = pu_vec[i] / wp;
+ wn_vec[i] = pu_vec[i] / 3.0e8;
+ wl_vec[i] = two_pi / wn_vec[i];
+ ev_vec[i] = pu_vec[i] * evc;
+ fprintf(
+ output,
+ "%5d%13.4lE%13.4lE%13.4lE%13.4lE%13.4lE\n",
+ (i + 1),
+ pu_vec[i],
+ wn_vec[i],
+ wl_vec[i],
+ ev_vec[i],
+ xi_vec[i]
+ );
+ }
+ break;
+ case 4:
+ fprintf(output, " JXI EVS WNS WLS PUS XIV\n");
+ for (int i = 0; i < number_of_scales; i++) {
+ ev_vec[i] = scale_vec[i];
+ pu_vec[i] = ev_vec[i] / evc;
+ xi_vec[i] = pu_vec[i] / wp;
+ wn_vec[i] = pu_vec[i] / 3.0e8;
+ wl_vec[i] = two_pi / wn_vec[i];
+ fprintf(
+ output,
+ "%5d%13.4lE%13.4lE%13.4lE%13.4lE%13.4lE\n",
+ (i + 1),
+ ev_vec[i],
+ wn_vec[i],
+ wl_vec[i],
+ pu_vec[i],
+ xi_vec[i]
+ );
+ }
+ break;
+ default:
+ throw UnrecognizedConfigurationException(
+ "Wonrg parameter set: unrecognized scale variable type "
+ + reference_variable_name
+ );
+ break;
+ }
+ } else { // idfc < 0, Dielectric functions are at XIP and XI is scale for dimensions
+ double pu, wn;
+ xi_vec = scale_vec;
+ pu = xip * wp;
+ wn = pu / 3.0e8;
+ fprintf(output, " XIP WN WL PU EV\n");
+ fprintf(output, " %13.4lE", xip);
+ fprintf(output, "%13.4lE", wn);
+ fprintf(output, "%13.4lE", two_pi / wn);
+ fprintf(output, "%13.4lE", pu);
+ fprintf(output, "%13.4lE\n", pu * evc);
+ fprintf(output, " SCALE FACTORS XI\n");
+ for (int i = 0; i < number_of_scales; i++)
+ fprintf(output, "%5d%13.4lE\n", (i + 1), xi_vec[i]);
+ }
+ if (idfc != 0) {
+ fprintf(output, " DIELECTRIC CONSTANTS\n");
+ for (int i473 = 1; i473 <= number_of_spheres; i473++) {
+ if (iog_vec[i473 - 1] != i473) continue;
+ ici = (nshl_vec[i473 - 1] + 1) / 2;
+ if (i473 == 1) ici += ies;
+ fprintf(output, " SPHERE N.%4d\n", i473);
+ for (int ic472 = 0; ic472 < ici; ic472++) {
+ double dc0_real = dc0_matrix[ic472][i473 - 1][0].real(), dc0_img = dc0_matrix[ic472][i473 - 1][0].imag();
+ fprintf(output, "%5d %12.4lE%12.4lE\n", (ic472 + 1), dc0_real, dc0_img);
+ }
+ }
+ } else {
+ fprintf(output, " DIELECTRIC FUNCTIONS\n");
+ for (int i478 = 1; i478 <= number_of_spheres; i478++) {
+ if (iog_vec[i478 - 1] != i478) continue;
+ ici = (nshl_vec[i478 - 1] + 1) / 2;
+ if (i478 == 1) ici += ies;
+ fprintf(output, " SPHERE N.%4d\n", i478);
+ for (int ic477 = 1; ic477 <= ici; ic477++) {
+ fprintf(output, " NONTRANSITION LAYER N.%2d, SCALE = %3s\n", ic477, reference_variable_name.c_str());
+ for (int jxi476 = 0; jxi476 < number_of_scales; jxi476++) {
+ double dc0_real = dc0_matrix[ic477 - 1][i478 - 1][jxi476].real();
+ double dc0_img = dc0_matrix[ic477 - 1][i478 - 1][jxi476].imag();
+ fprintf(output, "%5d%12.4lE%12.4lE\n", (jxi476 + 1), dc0_real, dc0_img);
}
- if (idfc != 0) {
- fprintf(output, " DIELECTRIC CONSTANTS\n");
- for (int i473 = 1; i473 <= number_of_spheres; i473++) {
- if (iog_vec[i473 - 1] != i473) continue;
- ici = (nshl_vec[i473 - 1] + 1) / 2;
- if (i473 == 1) ici += ies;
- fprintf(output, " SPHERE N.%4d\n", i473);
- for (int ic472 = 0; ic472 < ici; ic472++) {
- double dc0_real = dc0_matrix[ic472][i473 - 1][0].real(), dc0_img = dc0_matrix[ic472][i473 - 1][0].imag();
- fprintf(output, "%5d %12.4lE%12.4lE\n", (ic472 + 1), dc0_real, dc0_img);
- }
- }
- } else {
- fprintf(output, " DIELECTRIC FUNCTIONS\n");
- for (int i478 = 1; i478 <= number_of_spheres; i478++) {
- if (iog_vec[i478 - 1] != i478) continue;
- ici = (nshl_vec[i478 - 1] + 1) / 2;
- if (i478 == 1) ici += ies;
- fprintf(output, " SPHERE N.%4d\n", i478);
- for (int ic477 = 1; ic477 <= ici; ic477++) {
- fprintf(output, " NONTRANSITION LAYER N.%2d, SCALE = %3s\n", ic477, reference_variable_name.c_str());
- for (int jxi476 = 0; jxi476 < number_of_scales; jxi476++) {
- double dc0_real = dc0_matrix[ic477 - 1][i478 - 1][jxi476].real();
- double dc0_img = dc0_matrix[ic477 - 1][i478 - 1][jxi476].imag();
- fprintf(output, "%5d%12.4lE%12.4lE\n", (jxi476 + 1), dc0_real, dc0_img);
- }
- }
- }
- }
- fclose(output);
+ }
+ }
+ }
+ fclose(output);
}