From e3f610295af22d17d30cf1ea03e333d86e4f3610 Mon Sep 17 00:00:00 2001
From: Giovanni La Mura <giovanni.lamura@inaf.it>
Date: Thu, 28 Mar 2024 11:49:59 +0100
Subject: [PATCH] Provide method based access to the configuration parameters
of SPHERE
---
src/include/Configuration.h | 69 +++++++--
src/libnptm/Configuration.cpp | 2 +-
src/sphere/sphere.cpp | 280 +++++++++++++++++++++-------------
3 files changed, 237 insertions(+), 114 deletions(-)
diff --git a/src/include/Configuration.h b/src/include/Configuration.h
index dbebc711..f64330bf 100644
--- a/src/include/Configuration.h
+++ b/src/include/Configuration.h
@@ -55,7 +55,7 @@ class GeometryConfiguration {
//! Temporary work-around to allow cluster() and sphere() peeking in.
friend void cluster(std::string, std::string, std::string);
friend int cluster_jxi488_cycle(int, ScattererConfiguration *, GeometryConfiguration *, C1 *, C1_AddOns *, C2 *, C3 *, C4 *, C6 *, C9 *, std::FILE *, std::string, double *, double **, dcomplex **, double **, dcomplex **, double ****, double **, dcomplex **, int, int, int, int, int, int, int, double *, double *, double **, dcomplex **, double **, dcomplex **, double **, dcomplex **, double *, std::fstream &, double **, double **, double **, double **, double *, double *, double *, int, int, np_int, int, int, int, int, int, int, double, double, double, double, double, double, double, double, double, double, double, double, double, double, double, double, double, double *, double *, double *, double *, double *, double *, double *, double *, double *, double *, double &, double &, double &, double &, double &, double, double, int, dcomplex, double, double, np_int, dcomplex **, int, int, Logger *);
- friend void sphere(std::string, std::string, std::string);
+
protected:
//! \brief Number of spherical components.
int number_of_spheres;
@@ -231,7 +231,7 @@ class ScattererConfiguration {
//! Temporary work-around to allow cluster() and sphere() peeking in.
friend void cluster(std::string, std::string, std::string);
friend int cluster_jxi488_cycle(int, ScattererConfiguration *, GeometryConfiguration *, C1 *, C1_AddOns *, C2 *, C3 *, C4 *, C6 *, C9 *, std::FILE *, std::string, double *, double **, dcomplex **, double **, dcomplex **, double ****, double **, dcomplex **, int, int, int, int, int, int, int, double *, double *, double **, dcomplex **, double **, dcomplex **, double **, dcomplex **, double *, std::fstream &, double **, double **, double **, double **, double *, double *, double *, int, int, np_int, int, int, int, int, int, int, double, double, double, double, double, double, double, double, double, double, double, double, double, double, double, double, double, double *, double *, double *, double *, double *, double *, double *, double *, double *, double *, double &, double &, double &, double &, double &, double, double, int, dcomplex, double, double, np_int, dcomplex **, int, int, Logger *);
- friend void sphere(std::string, std::string, std::string);
+
protected:
//! \brief Matrix of dielectric parameters with size [NON_TRANS_LAYERS x N_SPHERES x N_SCALES].
dcomplex ***dc0_matrix;
@@ -407,16 +407,45 @@ public:
*/
static ScattererConfiguration* from_dedfb(std::string file_name);
- /*! \brief Get the ID of the configuration of a sphere by the sphere's index.
+ /*! \brief Get the dielectric constant of a material for a specific wavelength.
*
- * This is a specialized function to access the ID of the configuration group a
- * sphere belongs to, through the index of the sphere.
+ * Dielectric constants are stored in a 3D complex matrix, whose dimensions map
+ * to [NUMBER_OF_CONFIGURATIONS x NUMBER_OF_SPHERES x NUMBER_OF_SCALES]. This
+ * function extracts such values from the matrix through their indices.
*
- * \param index: `int` Index of the ID to be retrieved.
- * \return id: `int` The desired identifier.
+ * \param i: `int` Index of the configuration.
+ * \param j: `int` Index of the sphere.
+ * \param k: `int` Index of the current scale.
+ * \return radius: `dcomplex` The requested dielectric constant.
+ */
+ dcomplex get_dielectric_constant(int i, int j, int k) { return dc0_matrix[i][j][k]; }
+
+ /*! \brief Get the ID of a configuration from the index of the sphere.
+ *
+ * This is a specialized function to a configuration ID through the index of
+ * the sphere it applies to.
+ *
+ * \return ID: `int` ID of the configuration to be applied.
*/
int get_iog(int index) { return iog_vec[index]; }
+ /*! \brief Get the address of the configuration ID vector.
+ *
+ * This is a specialized function to access the configuration ID vector.
+ *
+ * \return ptr: `int *` Pointer to the configuration index vector.
+ */
+ int *get_iog_vec() { return iog_vec; }
+
+ /*! \brief Get the address of the layer number vector.
+ *
+ * This is a specialized function to access the vector of layer numbers for the
+ * different configurations.
+ *
+ * \return ptr: `int *` Pointer to the vector of layer numbers.
+ */
+ int *get_nshl() { return nshl_vec; }
+
/*! \brief Get the value of a parameter by name.
*
* The proper way to access read-only parameters from outside a class is to define
@@ -431,13 +460,33 @@ public:
*/
double get_param(std::string param_name);
- /*! \brief Get the ID of a sphere by its index.
+ /*! \brief Get the radius of a sphere by its index.
+ *
+ * This is a specialized function to get the radius of a sphere through its
+ * index.
+ *
+ * \param index: `int` Index of the ID to be retrieved.
+ * \return radius: `double` The requested sphere radius.
+ */
+ double get_radius(int index) { return radii_of_spheres[index]; }
+
+ /*! \brief Get the value of a scale by its index.
+ *
+ * This is a specialized function to access a scale (generally a wavelength),
+ * through its index.
+ *
+ * \param index: `int` Index of the scale to be retrieved.
+ * \return scale: `double` The desired scale.
+ */
+ double get_rcf(int row, int column) { return rcf[row][column]; }
+
+ /*! \brief Get the reference variable name.
*
* This is a specialized function to get the name of the reference variable as a
* string.
*
- * \param index: `int` Index of the ID to be retrieved.
- * \return id: `int` The desired identifier.
+ * \return name: `string` The name of the variable used to calculate wavelength /
+ * size scaling.
*/
std::string get_reference() { return reference_variable_name; }
diff --git a/src/libnptm/Configuration.cpp b/src/libnptm/Configuration.cpp
index 8cde1b7e..ac5d7531 100644
--- a/src/libnptm/Configuration.cpp
+++ b/src/libnptm/Configuration.cpp
@@ -233,7 +233,7 @@ GeometryConfiguration* GeometryConfiguration::from_legacy(string file_name) {
return conf;
}
-GeometryConfiguration::get_param(std::string param_name) {
+double GeometryConfiguration::get_param(std::string param_name) {
double value;
if (param_name.compare("number_of_spheres") == 0) value = (double)number_of_spheres;
else if (param_name.compare("nsph") == 0) value = (double)number_of_spheres;
diff --git a/src/sphere/sphere.cpp b/src/sphere/sphere.cpp
index 94860567..866368b5 100644
--- a/src/sphere/sphere.cpp
+++ b/src/sphere/sphere.cpp
@@ -17,6 +17,10 @@
#include "../include/errors.h"
#endif
+#ifndef INCLUDE_LOGGING_H_
+#include "../include/logging.h"
+#endif
+
#ifndef INCLUDE_CONFIGURATION_H_
#include "../include/Configuration.h"
#endif
@@ -42,15 +46,18 @@ using namespace std;
* \param output_path: `string` Directory to write the output files in.
*/
void sphere(string config_file, string data_file, string output_path) {
+ Logger *logger = new Logger(LOG_INFO);
dcomplex arg, s0, tfsas;
double th, ph;
- printf("INFO: making legacy configuration...\n");
+ logger->log("INFO: making legacy configuration...\n");
ScattererConfiguration *sconf = NULL;
try {
sconf = ScattererConfiguration::from_dedfb(config_file);
} catch(const OpenConfigurationFileException &ex) {
- printf("\nERROR: failed to open scatterer configuration file.\n");
- printf("FILE: %s\n", ex.what());
+ logger->err("\nERROR: failed to open scatterer configuration file.\n");
+ string message = ex.what();
+ logger->err("FILE: " + message + "\n");
+ delete logger;
exit(1);
}
sconf->write_formatted(output_path + "/c_OEDFB");
@@ -60,12 +67,16 @@ void sphere(string config_file, string data_file, string output_path) {
try {
gconf = GeometryConfiguration::from_legacy(data_file);
} catch(const OpenConfigurationFileException &ex) {
- printf("\nERROR: failed to open geometry configuration file.\n");
- printf("FILE: %s\n", ex.what());
+ logger->err("\nERROR: failed to open geometry configuration file.\n");
+ string message = ex.what();
+ logger->err("FILE: " + message + "\n");
if (sconf != NULL) delete sconf;
+ delete logger;
exit(1);
}
- if (sconf->number_of_spheres == gconf->number_of_spheres) {
+ int s_nsph = (int)sconf->get_param("nsph");
+ int nsph = (int)gconf->get_param("nsph");
+ if (s_nsph == nsph) {
int isq, ibf;
double *argi, *args, *gaps;
double cost, sint, cosp, sinp;
@@ -103,12 +114,31 @@ void sphere(string config_file, string data_file, string output_path) {
double frx = 0.0, fry = 0.0, frz = 0.0;
double cfmp, cfsp, sfmp, sfsp;
int jw;
- int nsph = gconf->number_of_spheres;
- C1 *c1 = new C1(nsph, gconf->l_max, sconf->nshl_vec, sconf->iog_vec);
+ int l_max = gconf->get_param("l_max");
+ int *iog_vec = sconf->get_iog_vec();
+ int *nshl_vec = sconf->get_nshl();
+ C1 *c1 = new C1(nsph, l_max, nshl_vec, iog_vec);
for (int i = 0; i < nsph; i++) {
- c1->ros[i] = sconf->radii_of_spheres[i];
+ c1->ros[i] = sconf->get_radius(i);
}
- C2 *c2 = new C2(nsph, 5, gconf->npnt, gconf->npntts);
+ int npnt = (int)gconf->get_param("npnt");
+ int npntts = (int)gconf->get_param("npntts");
+ int in_pol = (int)gconf->get_param("in_pol");
+ int meridional_type = (int)gconf->get_param("meridional_type");
+ int jwtm = (int)gconf->get_param("jwtm");
+ double in_theta_start = gconf->get_param("in_theta_start");
+ double in_theta_step = gconf->get_param("in_theta_step");
+ double in_theta_end = gconf->get_param("in_theta_end");
+ double sc_theta_start = gconf->get_param("sc_theta_start");
+ double sc_theta_step = gconf->get_param("sc_theta_step");
+ double sc_theta_end = gconf->get_param("sc_theta_end");
+ double in_phi_start = gconf->get_param("in_phi_start");
+ double in_phi_step = gconf->get_param("in_phi_step");
+ double in_phi_end = gconf->get_param("in_phi_end");
+ double sc_phi_start = gconf->get_param("sc_phi_start");
+ double sc_phi_step = gconf->get_param("sc_phi_step");
+ double sc_phi_end = gconf->get_param("sc_phi_end");
+ C2 *c2 = new C2(nsph, 5, npnt, npntts);
argi = new double[1];
args = new double[1];
gaps = new double[2];
@@ -117,41 +147,37 @@ void sphere(string config_file, string data_file, string output_path) {
fprintf(
output,
" %5d%5d%5d%5d%5d%5d\n",
- gconf->number_of_spheres,
- gconf->l_max,
- gconf->in_pol,
- gconf->npnt,
- gconf->npntts,
- gconf->meridional_type
+ nsph,
+ l_max,
+ in_pol,
+ npnt,
+ npntts,
+ meridional_type
);
fprintf(output, " READ(IR,*)TH,THSTP,THLST,THS,THSSTP,THSLST\n");
fprintf(
output,
" %9.3lE %9.3lE %9.3lE %9.3lE %9.3lE %9.3lE\n",
- gconf->in_theta_start,
- gconf->in_theta_step,
- gconf->in_theta_end,
- gconf->sc_theta_start,
- gconf->sc_theta_step,
- gconf->sc_theta_end
+ in_theta_start,
+ in_theta_step,
+ in_theta_end,
+ sc_theta_start,
+ sc_theta_step,
+ sc_theta_end
);
fprintf(output, " READ(IR,*)PH,PHSTP,PHLST,PHS,PHSSTP,PHSLST\n");
fprintf(
output,
" %9.3lE %9.3lE %9.3lE %9.3lE %9.3lE %9.3lE\n",
- gconf->in_phi_start,
- gconf->in_phi_step,
- gconf->in_phi_end,
- gconf->sc_phi_start,
- gconf->sc_phi_step,
- gconf->sc_phi_end
+ in_phi_start,
+ in_phi_step,
+ in_phi_end,
+ sc_phi_start,
+ sc_phi_step,
+ sc_phi_end
);
fprintf(output, " READ(IR,*)JWTM\n");
- fprintf(
- output,
- " %5d\n",
- gconf->jwtm
- );
+ fprintf(output, " %5d\n", jwtm);
fprintf(output, " READ(ITIN)NSPHT\n");
fprintf(output, " READ(ITIN)(IOG(I),I=1,NSPH)\n");
fprintf(output, " READ(ITIN)EXDC,WP,XIP,IDFC,NXI\n");
@@ -160,36 +186,36 @@ void sphere(string config_file, string data_file, string output_path) {
fprintf(output, " READ(ITIN)(RCF(I,NS),NS=1,NSH)\n \n");
double sml = 1.0e-3;
int nth = 0, nph = 0;
- if (gconf->in_theta_step != 0.0)
- nth = int((gconf->in_theta_end - gconf->in_theta_start) / gconf->in_theta_step + sml);
+ if (in_theta_step != 0.0)
+ nth = int((in_theta_end - in_theta_start) / in_theta_step + sml);
nth += 1;
- if (gconf->in_phi_step != 0.0)
- nph = int((gconf->in_phi_end - gconf->in_phi_start) / gconf->in_phi_step + sml);
+ if (in_phi_step != 0.0)
+ nph = int((in_phi_end - in_phi_start) / in_phi_step + sml);
nph += 1;
int nths = 0, nphs = 0;
double thsca;
- if (gconf->meridional_type > 1) { // ISAM > 1, fixed scattering angle
+ if (meridional_type > 1) { // ISAM > 1, fixed scattering angle
nths = 1;
- thsca = gconf->sc_theta_start - gconf->in_theta_start;
+ thsca = sc_theta_start - in_theta_start;
} else { //ISAM <= 1
- if (gconf->in_theta_step != 0.0)
- nths = int((gconf->sc_theta_end - gconf->sc_theta_start) / gconf->sc_theta_step + sml);
+ if (in_theta_step != 0.0)
+ nths = int((sc_theta_end - sc_theta_start) / sc_theta_step + sml);
nths += 1;
- if (gconf->meridional_type == 1) { // ISAM = 1
+ if (meridional_type == 1) { // ISAM = 1
nphs = 1;
} else { // ISAM < 1
- if (gconf->sc_phi_step != 0.0)
- nphs = int((gconf->sc_phi_end - gconf->sc_phi_start) / gconf->sc_phi_step + sml);
+ if (sc_phi_step != 0.0)
+ nphs = int((sc_phi_end - sc_phi_start) / sc_phi_step + sml);
nphs += 1;
}
}
int nk = nth * nph;
int nks = nths * nphs;
int nkks = nk * nks;
- double th1 = gconf->in_theta_start;
- double ph1 = gconf->in_phi_start;
- double ths1 = gconf->sc_theta_start;
- double phs1 = gconf->sc_phi_start;
+ double th1 = in_theta_start;
+ double ph1 = in_phi_start;
+ double ths1 = sc_theta_start;
+ double phs1 = sc_phi_start;
const double half_pi = acos(0.0);
const double pi = 2.0 * half_pi;
double gcs = 0.0;
@@ -201,13 +227,13 @@ void sphere(string config_file, string data_file, string output_path) {
c1->gcsv[i116] = gcss;
int nsh = c1->nshl[i116];
for (int j115 = 0; j115 < nsh; j115++) {
- c1->rc[i116][j115] = sconf->rcf[i116][j115] * c1->ros[i116];
+ c1->rc[i116][j115] = sconf->get_rcf(i116, j115) * c1->ros[i116];
}
}
gcs += c1->gcsv[iogi];
}
- double ****zpv = new double***[gconf->l_max]; //[gconf->l_max][3][2][2]; // Matrix: dim[LM x 3 x 2 x 2]
- for (int zi = 0; zi < gconf->l_max; zi++) {
+ double ****zpv = new double***[l_max]; //[l_max][3][2][2]; // Matrix: dim[LM x 3 x 2 x 2]
+ for (int zi = 0; zi < l_max; zi++) {
zpv[zi] = new double**[3];
for (int zj = 0; zj < 3; zj++) {
zpv[zi][zj] = new double*[2];
@@ -216,8 +242,9 @@ void sphere(string config_file, string data_file, string output_path) {
}
}
}
- thdps(gconf->l_max, zpv);
- double exri = sqrt(sconf->exdc);
+ thdps(l_max, zpv);
+ double exdc = sconf->get_param("exdc");
+ double exri = sqrt(exdc);
fprintf(output, " REFR. INDEX OF EXTERNAL MEDIUM=%15.7lE\n", exri);
fstream tppoan;
string tppoan_name = output_path + "/c_TPPOAN";
@@ -225,9 +252,9 @@ void sphere(string config_file, string data_file, string output_path) {
if (tppoan.is_open()) {
int imode = 10;
tppoan.write(reinterpret_cast<char *>(&imode), sizeof(int));
- tppoan.write(reinterpret_cast<char *>(&(gconf->meridional_type)), sizeof(int));
- tppoan.write(reinterpret_cast<char *>(&(gconf->in_pol)), sizeof(int));
- tppoan.write(reinterpret_cast<char *>(&(sconf->number_of_scales)), sizeof(int));
+ tppoan.write(reinterpret_cast<char *>(&(meridional_type)), sizeof(int));
+ tppoan.write(reinterpret_cast<char *>(&(in_pol)), sizeof(int));
+ tppoan.write(reinterpret_cast<char *>(&s_nsph), sizeof(int));
tppoan.write(reinterpret_cast<char *>(&nth), sizeof(int));
tppoan.write(reinterpret_cast<char *>(&nph), sizeof(int));
tppoan.write(reinterpret_cast<char *>(&nths), sizeof(int));
@@ -235,24 +262,28 @@ void sphere(string config_file, string data_file, string output_path) {
tppoan.write(reinterpret_cast<char *>(&nsph), sizeof(int));
for (int nsi = 0; nsi < nsph; nsi++)
- tppoan.write(reinterpret_cast<char *>(&(sconf->iog_vec[nsi])), sizeof(int));
- if (gconf->in_pol == 0) fprintf(output, " LIN\n");
+ tppoan.write(reinterpret_cast<char *>(&(iog_vec[nsi])), sizeof(int));
+ if (in_pol == 0) fprintf(output, " LIN\n");
else fprintf(output, " CIRC\n");
fprintf(output, " \n");
- double wn = sconf->wp / 3.0e8;
+ double wp = sconf->get_param("wp");
+ double xip = sconf->get_param("xip");
+ double wn = wp / 3.0e8;
double sqsfi = 1.0;
double vk, vkarg;
- if (sconf->idfc < 0) {
- vk = sconf->xip * wn;
+ int idfc = (int)sconf->get_param("idfc");
+ int nxi = (int)sconf->get_param("number_of_scales");
+ if (idfc < 0) {
+ vk = xip * wn;
fprintf(output, " VK=%15.7lE, XI IS SCALE FACTOR FOR LENGTHS\n", vk);
fprintf(output, " \n");
}
- for (int jxi488 = 0; jxi488 < sconf->number_of_scales; jxi488++) {
+ for (int jxi488 = 0; jxi488 < nxi; jxi488++) {
int jxi = jxi488 + 1;
- printf("INFO: running scale iteration %d of %d...\n", jxi, sconf->number_of_scales);
+ logger->log("INFO: running scale iteration " + to_string(jxi) + " of " + to_string(nxi) + ".\n");
fprintf(output, "========== JXI =%3d ====================\n", jxi);
- double xi = sconf->scale_vec[jxi488];
- if (sconf->idfc >= 0) {
+ double xi = sconf->get_scale(jxi488);
+ if (idfc >= 0) {
vk = xi * wn;
vkarg = vk;
fprintf(output, " VK=%15.7lE, XI=%15.7lE\n", xi, vk);
@@ -264,54 +295,96 @@ void sphere(string config_file, string data_file, string output_path) {
tppoan.write(reinterpret_cast<char *>(&vk), sizeof(double));
for (int i132 = 0; i132 < nsph; i132++) {
int i = i132 + 1;
- int iogi = sconf->iog_vec[i132];
+ int iogi = iog_vec[i132];
if (iogi != i) {
- for (int l123 = 0; l123 < gconf->l_max; l123++) {
+ for (int l123 = 0; l123 < l_max; l123++) {
c1->rmi[l123][i132] = c1->rmi[l123][iogi - 1];
c1->rei[l123][i132] = c1->rei[l123][iogi - 1];
}
continue; // i132
}
- int nsh = sconf->nshl_vec[i132];
+ int nsh = nshl_vec[i132];
int ici = (nsh + 1) / 2;
- if (sconf->idfc == 0) {
+ if (idfc == 0) {
for (int ic = 0; ic < ici; ic++)
- c2->dc0[ic] = sconf->dc0_matrix[ic][i132][0]; // WARNING: IDFC=0 is not tested!
+ c2->dc0[ic] = sconf->get_dielectric_constant(ic, i132, 0); // WARNING: IDFC=0 is not tested!
} else { // IDFC != 0
if (jxi == 1) {
for (int ic = 0; ic < ici; ic++) {
- c2->dc0[ic] = sconf->dc0_matrix[ic][i132][jxi488];
+ c2->dc0[ic] = sconf->get_dielectric_constant(ic, i132, jxi488);
}
}
}
- if (nsh % 2 == 0) c2->dc0[ici] = sconf->exdc;
+ if (nsh % 2 == 0) c2->dc0[ici] = exdc;
int jer = 0;
int lcalc = 0;
- dme(
- gconf->l_max, i, gconf->npnt, gconf->npntts, vkarg,
- sconf->exdc, exri, c1, c2, jer, lcalc, arg
- );
+ dme(l_max, i, npnt, npntts, vkarg, exdc, exri, c1, c2, jer, lcalc, arg);
if (jer != 0) {
fprintf(output, " STOP IN DME\n");
fprintf(output, " AT %1d LCALC=%3d TOO SMALL WITH ARG=%15.7lE+i(%15.7lE)\n", jer, lcalc, real(arg), imag(arg));
tppoan.close();
fclose(output);
+ delete sconf;
+ delete gconf;
+ delete c1;
+ delete c2;
+ for (int zi = l_max - 1; zi > -1; zi--) {
+ for (int zj = 0; zj < 3; zj++) {
+ for (int zk = 0; zk < 2; zk++) {
+ delete[] zpv[zi][zj][zk];
+ }
+ delete[] zpv[zi][zj];
+ }
+ delete[] zpv[zi];
+ }
+ delete[] zpv;
+ delete[] duk;
+ delete[] u;
+ delete[] us;
+ delete[] un;
+ delete[] uns;
+ delete[] up;
+ delete[] ups;
+ delete[] upmp;
+ delete[] upsmp;
+ delete[] unmp;
+ delete[] unsmp;
+ delete[] argi;
+ delete[] args;
+ delete[] gaps;
+ for (int i = 3; i > -1; i--) {
+ delete[] cmul[i];
+ delete[] cmullr[i];
+ }
+ delete[] cmul;
+ delete[] cmullr;
+ for (int ti = 1; ti > -1; ti--) {
+ delete[] tqse[ti];
+ delete[] tqss[ti];
+ delete[] tqspe[ti];
+ delete[] tqsps[ti];
+ }
+ delete[] tqse;
+ delete[] tqss;
+ delete[] tqspe;
+ delete[] tqsps;
+ delete logger;
return;
}
} // i132
- if (sconf->idfc >= 0 and nsph == 1 and jxi == gconf->jwtm) {
+ if (idfc >= 0 and nsph == 1 and jxi == jwtm) {
// This is the condition that writes the transition matrix to output.
- TransitionMatrix ttms(gconf->l_max, vk, exri, c1->rmi, c1->rei, sconf->radii_of_spheres[0]);
+ TransitionMatrix ttms(l_max, vk, exri, c1->rmi, c1->rei, sconf->get_radius(0));
string ttms_name = output_path + "/c_TTMS.hd5";
ttms.write_binary(ttms_name, "HDF5");
ttms_name = output_path + "/c_TTMS";
ttms.write_binary(ttms_name);
}
double cs0 = 0.25 * vk * vk * vk / half_pi;
- sscr0(tfsas, nsph, gconf->l_max, vk, exri, c1);
- double sqk = vk * vk * sconf->exdc;
- aps(zpv, gconf->l_max, nsph, c1, sqk, gaps);
- rabas(gconf->in_pol, gconf->l_max, nsph, c1, tqse, tqspe, tqss, tqsps);
+ sscr0(tfsas, nsph, l_max, vk, exri, c1);
+ double sqk = vk * vk * exdc;
+ aps(zpv, l_max, nsph, c1, sqk, gaps);
+ rabas(in_pol, l_max, nsph, c1, tqse, tqspe, tqss, tqsps);
for (int i170 = 0; i170 < nsph; i170++) {
int i = i170 + 1;
if (c1->iog[i170] >= i) {
@@ -406,8 +479,8 @@ void sphere(string config_file, string data_file, string output_path) {
if (!goto182) {
upvmp(th, ph, 0, cost, sint, cosp, sinp, u, upmp, unmp);
}
- if (gconf->meridional_type >= 0) {
- wmamp(0, cost, sint, cosp, sinp, gconf->in_pol, gconf->l_max, 0, nsph, argi, u, upmp, unmp, c1);
+ if (meridional_type >= 0) {
+ wmamp(0, cost, sint, cosp, sinp, in_pol, l_max, 0, nsph, argi, u, upmp, unmp, c1);
for (int i183 = 0; i183 < nsph; i183++) {
double rapr = c1->sexs[i183] - gaps[i183];
frx = rapr * u[0];
@@ -422,8 +495,8 @@ void sphere(string config_file, string data_file, string output_path) {
int jths = jths482 + 1;
double ths = thsl;
int icspnv = 0;
- if (gconf->meridional_type > 1) ths = th + thsca;
- if (gconf->meridional_type >= 1) {
+ if (meridional_type > 1) ths = th + thsca;
+ if (meridional_type >= 1) {
phsph = 0.0;
if ((ths < 0.0) || (ths > 180.0)) phsph = 180.0;
if (ths < 0.0) ths *= -1.0;
@@ -433,24 +506,24 @@ void sphere(string config_file, string data_file, string output_path) {
double phs = phs1;
for (int jphs480 = 0; jphs480 < nphs; jphs480++) {
int jphs = jphs480 + 1;
- if (gconf->meridional_type >= 1) {
+ if (meridional_type >= 1) {
phs = ph + phsph;
if (phs >= 360.0) phs -= 360.0;
}
bool goto190 = (nks == 1) && ((jxi > 1) || (jth > 1) || (jph > 1));
if (!goto190) {
upvmp(ths, phs, icspnv, costs, sints, cosps, sinps, us, upsmp, unsmp);
- if (gconf->meridional_type >= 0) {
- wmamp(2, costs, sints, cosps, sinps, gconf->in_pol, gconf->l_max, 0, nsph, args, us, upsmp, unsmp, c1);
+ if (meridional_type >= 0) {
+ wmamp(2, costs, sints, cosps, sinps, in_pol, l_max, 0, nsph, args, us, upsmp, unsmp, c1);
}
}
if (nkks != 0 || jxi == 1) {
upvsp(u, upmp, unmp, us, upsmp, unsmp, up, un, ups, uns, duk, isq, ibf, scan, cfmp, sfmp, cfsp, sfsp);
- if (gconf->meridional_type < 0) {
+ if (meridional_type < 0) {
wmasp(
cost, sint, cosp, sinp, costs, sints, cosps, sinps,
- u, up, un, us, ups, uns, isq, ibf, gconf->in_pol,
- gconf->l_max, 0, nsph, argi, args, c1
+ u, up, un, us, ups, uns, isq, ibf, in_pol,
+ l_max, 0, nsph, argi, args, c1
);
}
for (int i193 = 0; i193 < 3; i193++) {
@@ -458,7 +531,7 @@ void sphere(string config_file, string data_file, string output_path) {
uns[i193] = unsmp[i193];
}
}
- if (gconf->meridional_type < 0) jw = 1;
+ if (meridional_type < 0) jw = 1;
tppoan.write(reinterpret_cast<char *>(&th), sizeof(double));
tppoan.write(reinterpret_cast<char *>(&ph), sizeof(double));
tppoan.write(reinterpret_cast<char *>(&ths), sizeof(double));
@@ -483,13 +556,13 @@ void sphere(string config_file, string data_file, string output_path) {
fprintf(output, " SCAND=%10.3lE\n", scan);
fprintf(output, " CFMP=%15.7lE, SFMP=%15.7lE\n", cfmp, sfmp);
fprintf(output, " CFSP=%15.7lE, SFSP=%15.7lE\n", cfsp, sfsp);
- if (gconf->meridional_type >= 0) {
+ if (meridional_type >= 0) {
fprintf(output, " UNI=(%12.5lE,%12.5lE,%12.5lE)\n", un[0], un[1], un[2]);
fprintf(output, " UNS=(%12.5lE,%12.5lE,%12.5lE)\n", uns[0], uns[1], uns[2]);
} else {
fprintf(output, " UN=(%12.5lE,%12.5lE,%12.5lE)\n", un[0], un[1], un[2]);
}
- sscr2(nsph, gconf->l_max, vk, exri, c1);
+ sscr2(nsph, l_max, vk, exri, c1);
for (int ns226 = 0; ns226 < nsph; ns226++) {
int ns = ns226 + 1;
fprintf(output, " SPHERE %2d\n", ns);
@@ -538,24 +611,24 @@ void sphere(string config_file, string data_file, string output_path) {
}
}
} // ns226 loop
- if (gconf->meridional_type < 1) phs += gconf->sc_phi_step;
+ if (meridional_type < 1) phs += sc_phi_step;
} // jphs480 loop
- if (gconf->meridional_type <= 1) thsl += gconf->sc_theta_step;
+ if (meridional_type <= 1) thsl += sc_theta_step;
} // jths482 loop
- ph += gconf->in_phi_step;
+ ph += in_phi_step;
} // jph484 loop on elevation
- th += gconf->in_theta_step;
+ th += in_theta_step;
} // jth486 loop on azimuth
- printf("INFO: done scale.\n");
+ logger->log("INFO: finished scale iteration " + to_string(jxi) + " of " + to_string(nxi) + ".\n");
} //jxi488 loop on scales
tppoan.close();
} else { // In case TPPOAN could not be opened. Should never happen.
- printf("ERROR: failed to open TPPOAN file.\n");
+ logger->err("ERROR: failed to open TPPOAN file.\n");
}
fclose(output);
delete c1;
delete c2;
- for (int zi = gconf->l_max - 1; zi > -1; zi--) {
+ for (int zi = l_max - 1; zi > -1; zi--) {
for (int zj = 0; zj < 3; zj++) {
for (int zk = 0; zk < 2; zk++) {
delete[] zpv[zi][zj][zk];
@@ -595,7 +668,7 @@ void sphere(string config_file, string data_file, string output_path) {
delete[] tqss;
delete[] tqspe;
delete[] tqsps;
- printf("Finished: output written to %s.\n", (output_path + "/c_OSPH").c_str());
+ logger->log("Finished: output written to " + output_path + "/c_OSPH.\n");
} else { // NSPH mismatch between geometry and scatterer configurations.
throw UnrecognizedConfigurationException(
"Inconsistent geometry and scatterer configurations."
@@ -603,4 +676,5 @@ void sphere(string config_file, string data_file, string output_path) {
}
delete sconf;
delete gconf;
+ delete logger;
}
--
GitLab