From a6765e93b96a3fa83c57274512be1b7e623ba107 Mon Sep 17 00:00:00 2001
From: Randolf Beerwerth <r.beerwerth@fz-juelich.de>
Date: Wed, 16 Sep 2020 14:43:16 +0200
Subject: [PATCH] make types more specific, some cleanup and enhancements

---
 .../SpecularMagneticNewNCStrategy.cpp         | 40 ++++++-------
 .../SpecularMagneticNewStrategy.cpp           |  8 +--
 .../SpecularMagneticNewTanhStrategy.cpp       | 33 +++++-----
 Core/RT/MatrixRTCoefficients_v3.cpp           | 60 +++++++------------
 4 files changed, 59 insertions(+), 82 deletions(-)

diff --git a/Core/Multilayer/SpecularMagneticNewNCStrategy.cpp b/Core/Multilayer/SpecularMagneticNewNCStrategy.cpp
index da03ec2b0ec..426fc0c91ff 100644
--- a/Core/Multilayer/SpecularMagneticNewNCStrategy.cpp
+++ b/Core/Multilayer/SpecularMagneticNewNCStrategy.cpp
@@ -24,51 +24,49 @@ SpecularMagneticNewNCStrategy::computeRoughnessMatrices(const MatrixRTCoefficien
                                                         const MatrixRTCoefficients_v3& coeff_i1,
                                                         double sigma) const
 {
-    auto beta_i  = coeff_i.m_lambda(1) - coeff_i.m_lambda(0);
-    auto beta_i1 = coeff_i1.m_lambda(1) - coeff_i1.m_lambda(0);
+    complex_t beta_i  = coeff_i.m_lambda(1) - coeff_i.m_lambda(0);
+    complex_t beta_i1 = coeff_i1.m_lambda(1) - coeff_i1.m_lambda(0);
 
     auto roughness_matrix = [sigma, &coeff_i, &coeff_i1, beta_i, beta_i1](double sign){
 
-        auto alpha_p = coeff_i1.m_lambda(0) + coeff_i1.m_lambda(1) +
+        const complex_t alpha_p = coeff_i1.m_lambda(0) + coeff_i1.m_lambda(1) +
                         sign * (coeff_i.m_lambda(0) + coeff_i.m_lambda(1));
         auto b_p_vec = beta_i1 * coeff_i1.m_b + sign * beta_i * coeff_i.m_b;
 
         auto square = [](auto & v){
             return v.x() * v.x() + v.y() * v.y() + v.z() * v.z();
         };
-        auto beta_p = std::sqrt(checkForUnderflow(square(b_p_vec)));
+        complex_t beta_p = std::sqrt(checkForUnderflow(square(b_p_vec)));
         b_p_vec /= beta_p;
 
-        auto alpha_pp = -(alpha_p * alpha_p + beta_p * beta_p) * sigma * sigma / 8.;
-        auto beta_pp  = - alpha_p * beta_p * sigma * sigma / 4.;
+        const complex_t alpha_pp = -(alpha_p * alpha_p + beta_p * beta_p) * sigma * sigma / 8.;
+        const complex_t beta_pp  = - alpha_p * beta_p * sigma * sigma / 4.;
 
         Eigen::Matrix2cd QL, QR;
 
-        auto factor1 = std::sqrt(2. * (1. + b_p_vec.z()));
-        auto factor2 = std::sqrt(2. * (1. - b_p_vec.z()));
+        const complex_t factor1 = std::sqrt(2. * (1. + b_p_vec.z()));
+        const complex_t factor2 = std::sqrt(2. * (1. - b_p_vec.z()));
         QL << (b_p_vec.z() + 1.)/factor1,                (b_p_vec.z() - 1.)/factor2,
                 (b_p_vec.x() + I * b_p_vec.y())/factor1, (b_p_vec.x() + I * b_p_vec.y())/factor2;
         QR << (b_p_vec.z() + 1.)/factor1,   (b_p_vec.x() - I * b_p_vec.y())/factor1,
                 (b_p_vec.z() - 1.)/factor2, (b_p_vec.x() - I * b_p_vec.y())/factor2;
 
-        auto exp1 = Eigen::Matrix2cd(Eigen::DiagonalMatrix<complex_t, 2>(
-                            {std::exp(alpha_pp), std::exp(alpha_pp)}));
+        const Eigen::Matrix2cd exp1 = Eigen::DiagonalMatrix<complex_t, 2>(
+                            {std::exp(alpha_pp), std::exp(alpha_pp)});
 
-        Eigen::Matrix2cd roughnessMatrix;
         if( std::abs(beta_p) > std::numeric_limits<double>::epsilon() * 10. )
         {
             Eigen::Matrix2cd exp2 = Eigen::Matrix2cd(Eigen::DiagonalMatrix<complex_t, 2>(
                                 {std::exp(beta_pp), std::exp(-beta_pp)}));
-            roughnessMatrix = exp1 * QL * exp2 * QR;
+            return Eigen::Matrix2cd{exp1 * QL * exp2 * QR};
+        }
 
-        }else
-            roughnessMatrix = exp1;
 
-        return roughnessMatrix;
+        return exp1;
     };
 
-    auto roughness_sum  = roughness_matrix(1.);
-    auto roughness_diff = roughness_matrix(-1.);
+    const Eigen::Matrix2cd roughness_sum  = roughness_matrix(1.);
+    const Eigen::Matrix2cd roughness_diff = roughness_matrix(-1.);
 
     return {roughness_sum, roughness_diff};
 }
@@ -81,14 +79,14 @@ SpecularMagneticNewNCStrategy::computeBackwardsSubmatrices(const MatrixRTCoeffic
     Eigen::Matrix2cd roughness_sum{Eigen::Matrix2cd::Identity()};
     Eigen::Matrix2cd roughness_diff{Eigen::Matrix2cd::Identity()};
     if (sigma != 0.) {
-        auto ret = computeRoughnessMatrices(coeff_i, coeff_i1, sigma);
+        const auto ret = computeRoughnessMatrices(coeff_i, coeff_i1, sigma);
         roughness_sum = std::get<0>(ret);
         roughness_diff = std::get<1>(ret);
     }
 
-    auto P = Eigen::Matrix2cd(coeff_i.computeInverseP() * coeff_i1.computeP());
-    auto mp = 0.5 * (Eigen::Matrix2cd::Identity() + P) * roughness_diff;
-    auto mm = 0.5 * (Eigen::Matrix2cd::Identity() - P) * roughness_sum;
+    const auto P = Eigen::Matrix2cd(coeff_i.computeInverseP() * coeff_i1.computeP());
+    const auto mp = 0.5 * (Eigen::Matrix2cd::Identity() + P) * roughness_diff;
+    const auto mm = 0.5 * (Eigen::Matrix2cd::Identity() - P) * roughness_sum;
 
     return {mp, mm};
 }
diff --git a/Core/Multilayer/SpecularMagneticNewStrategy.cpp b/Core/Multilayer/SpecularMagneticNewStrategy.cpp
index 53d41f45cf0..5462c664390 100644
--- a/Core/Multilayer/SpecularMagneticNewStrategy.cpp
+++ b/Core/Multilayer/SpecularMagneticNewStrategy.cpp
@@ -97,7 +97,7 @@ SpecularMagneticNewStrategy::computeTR(const std::vector<Slice>& slices,
 void SpecularMagneticNewStrategy::calculateUpwards(
     std::vector<MatrixRTCoefficients_v3>& coeff, const std::vector<Slice>& slices) const
 {
-    auto N = slices.size();
+    const auto N = slices.size();
     std::vector<Eigen::Matrix2cd> SMatrices(N-1);
     std::vector<complex_t> Normalization(N-1);
 
@@ -111,7 +111,7 @@ void SpecularMagneticNewStrategy::calculateUpwards(
             sigma = roughness->getSigma();
 
         // compute the 2x2 submatrices in the write-up denoted as P+, P- and delta
-        auto mpmm = computeBackwardsSubmatrices(coeff[i], coeff[i + 1], sigma);
+        const auto mpmm = computeBackwardsSubmatrices(coeff[i], coeff[i + 1], sigma);
         const Eigen::Matrix2cd mp = mpmm.first;
         const Eigen::Matrix2cd mm = mpmm.second;
 
@@ -122,7 +122,7 @@ void SpecularMagneticNewStrategy::calculateUpwards(
         Si = mp + mm * coeff[i+1].m_R;
         S << Si(1, 1) , -Si(0, 1),
              -Si(1, 0), Si(0, 0);
-        auto norm = S(0, 0) * S(1, 1) - S(0, 1) * S(1, 0);
+        const complex_t norm = S(0, 0) * S(1, 1) - S(0, 1) * S(1, 0);
         S = S * delta;
 
         // store the rotation matrix and normalization constant in order to rotate
@@ -142,7 +142,7 @@ void SpecularMagneticNewStrategy::calculateUpwards(
     // normalization constants. In addition rotate the polarization by the amount
     // that was rotated above the current interface
     // if the normalization overflows, all amplitudes below that point are set to zero
-    auto dumpingFactor = complex_t(1, 0);
+    complex_t dumpingFactor = 1;
     Eigen::Matrix2cd S = Eigen::Matrix2cd::Identity();
     for (size_t i = 1; i < N; ++i) {
         dumpingFactor = dumpingFactor * Normalization[i - 1];
diff --git a/Core/Multilayer/SpecularMagneticNewTanhStrategy.cpp b/Core/Multilayer/SpecularMagneticNewTanhStrategy.cpp
index 4725681b332..7ca12455e59 100644
--- a/Core/Multilayer/SpecularMagneticNewTanhStrategy.cpp
+++ b/Core/Multilayer/SpecularMagneticNewTanhStrategy.cpp
@@ -28,41 +28,38 @@ SpecularMagneticNewTanhStrategy::computeRoughnessMatrix(const MatrixRTCoefficien
     if (sigma < 10 * std::numeric_limits<double>::epsilon())
         return Eigen::Matrix2cd{Eigen::Matrix2cd::Identity()};
 
-    const auto sigeff = pi2_15 * sigma;
-    auto b = coeff.m_b;
+    const double sigeff = pi2_15 * sigma;
+    const auto b = coeff.m_b;
 
-    Eigen::Matrix2cd R;
     if (std::abs(b.mag() - 1.) < std::numeric_limits<double>::epsilon() * 10.) {
         Eigen::Matrix2cd Q;
-        auto factor1 = 2. * (1. + b.z());
+        const double factor1 = 2. * (1. + b.z());
         Q << (1. + b.z()),    (I * b.y() - b.x()),
                 (b.x() + I * b.y()), (b.z() + 1.);
 
-        auto l1 = std::sqrt(MathFunctions::tanhc(sigeff * coeff.m_lambda(1)));
-        auto l2 = std::sqrt(MathFunctions::tanhc(sigeff * coeff.m_lambda(0)));
+        complex_t l1 = std::sqrt(MathFunctions::tanhc(sigeff * coeff.m_lambda(1)));
+        complex_t l2 = std::sqrt(MathFunctions::tanhc(sigeff * coeff.m_lambda(0)));
 
         if (inverse) {
             l1 = 1. / l1;
             l2 = 1. / l2;
         }
 
-        auto lambda = Eigen::Matrix2cd(Eigen::DiagonalMatrix<complex_t, 2>({l1, l2}));
+        const Eigen::Matrix2cd lambda = Eigen::DiagonalMatrix<complex_t, 2>({l1, l2});
 
-        R = Q * lambda * Q.adjoint() / factor1;
+        return Q * lambda * Q.adjoint() / factor1;
 
     } else if (b.mag() < 10 * std::numeric_limits<double>::epsilon()) {
-        auto alpha =
+        complex_t alpha =
             std::sqrt(MathFunctions::tanhc(0.5 * sigeff * (coeff.m_lambda(1) + coeff.m_lambda(0))));
         if (inverse)
             alpha = 1. / alpha;
-        auto lambda = Eigen::Matrix2cd(Eigen::DiagonalMatrix<complex_t, 2>({alpha, alpha}));
+        const Eigen::Matrix2cd lambda = Eigen::DiagonalMatrix<complex_t, 2>({alpha, alpha});
 
-        R = lambda;
-
-    } else
-        throw std::runtime_error("Broken magnetic field vector");
+        return lambda;
+    }
 
-    return R;
+    throw std::runtime_error("Broken magnetic field vector");
 }
 
 std::pair<Eigen::Matrix2cd, Eigen::Matrix2cd>
@@ -81,9 +78,9 @@ SpecularMagneticNewTanhStrategy::computeBackwardsSubmatrices(
                * computeRoughnessMatrix(coeff_i1, sigma, true);
     }
 
-    auto mproduct = Eigen::Matrix2cd(coeff_i.computeInverseP() * coeff_i1.computeP());
-    auto mp = 0.5 * (RInv + mproduct * R);
-    auto mm = 0.5 * (RInv - mproduct * R);
+    const Eigen::Matrix2cd mproduct = coeff_i.computeInverseP() * coeff_i1.computeP();
+    const Eigen::Matrix2cd mp = 0.5 * (RInv + mproduct * R);
+    const Eigen::Matrix2cd mm = 0.5 * (RInv - mproduct * R);
 
     return {mp, mm};
 }
diff --git a/Core/RT/MatrixRTCoefficients_v3.cpp b/Core/RT/MatrixRTCoefficients_v3.cpp
index 32e2e74b8c8..18b4ec09c8c 100644
--- a/Core/RT/MatrixRTCoefficients_v3.cpp
+++ b/Core/RT/MatrixRTCoefficients_v3.cpp
@@ -14,7 +14,6 @@
 
 #include "Core/RT/MatrixRTCoefficients_v3.h"
 #include "Core/Basics/Assert.h"
-#include <iostream>
 
 namespace
 {
@@ -45,25 +44,21 @@ MatrixRTCoefficients_v3* MatrixRTCoefficients_v3::clone() const
 Eigen::Matrix2cd MatrixRTCoefficients_v3::TransformationMatrix(complex_t eigenvalue,
                                                                Eigen::Vector2d selection) const
 {
-    Eigen::Matrix2cd result;
-
-    auto exp2 = Eigen::Matrix2cd(Eigen::DiagonalMatrix<complex_t, 2>(selection));
+    const Eigen::Matrix2cd exp2 = Eigen::DiagonalMatrix<complex_t, 2>(selection);
 
     if (std::abs(m_b.mag() - 1.) < eps) {
         Eigen::Matrix2cd Q;
-        auto factor1 = 2. * (1. + m_b.z());
+        const double factor1 = 2. * (1. + m_b.z());
         Q << (1. + m_b.z()), (I * m_b.y() - m_b.x()), (m_b.x() + I * m_b.y()), (m_b.z() + 1.);
-        result = Q * exp2 * Q.adjoint() / factor1;
+        return Q * exp2 * Q.adjoint() / factor1;
 
     } else if (m_b.mag() < eps && eigenvalue != 0.)
-        result = exp2;
+        return exp2;
     //        result = Eigen::Matrix2cd(Eigen::DiagonalMatrix<complex_t, 2>({0.5, 0.5}));
     else if (m_b.mag() < eps && eigenvalue == 0.)
-        result = Eigen::Matrix2cd(Eigen::DiagonalMatrix<complex_t, 2>({0.5, 0.5}));
-    else
-        throw std::runtime_error("Broken magnetic field vector");
+        return Eigen::Matrix2cd(Eigen::DiagonalMatrix<complex_t, 2>({0.5, 0.5}));
 
-    return result;
+    throw std::runtime_error("Broken magnetic field vector");
 }
 
 Eigen::Matrix2cd MatrixRTCoefficients_v3::T1Matrix() const
@@ -78,52 +73,42 @@ Eigen::Matrix2cd MatrixRTCoefficients_v3::T2Matrix() const
 
 Eigen::Vector2cd MatrixRTCoefficients_v3::T1plus() const
 {
-    auto mat = T1Matrix();
-    return mat * m_T.col(0);
-    ;
+    return T1Matrix() * m_T.col(0);
 }
 
 Eigen::Vector2cd MatrixRTCoefficients_v3::R1plus() const
 {
-    auto mat = T1Matrix();
-    return mat * m_R.col(0);
-    ;
+    return T1Matrix() * m_R.col(0);
 }
 
 Eigen::Vector2cd MatrixRTCoefficients_v3::T2plus() const
 {
-    auto mat = T2Matrix();
-    return mat * m_T.col(0);
+    return T2Matrix() * m_T.col(0);
 }
 
 Eigen::Vector2cd MatrixRTCoefficients_v3::R2plus() const
 {
-    auto mat = T2Matrix();
-    return mat * m_R.col(0);
+    return T2Matrix() * m_R.col(0);
 }
 
 Eigen::Vector2cd MatrixRTCoefficients_v3::T1min() const
 {
-    auto mat = T1Matrix();
-    return mat * m_T.col(1);
+    return T1Matrix() * m_T.col(1);
 }
 
 Eigen::Vector2cd MatrixRTCoefficients_v3::R1min() const
 {
-    auto mat = T1Matrix();
-    return mat * m_R.col(1);
+    return T1Matrix() * m_R.col(1);
 }
 
 Eigen::Vector2cd MatrixRTCoefficients_v3::T2min() const
 {
-    auto mat = T2Matrix();
-    return mat * m_T.col(1);
+    return T2Matrix() * m_T.col(1);
 }
 
 Eigen::Vector2cd MatrixRTCoefficients_v3::R2min() const
 {
-    auto mat = T2Matrix();
-    return mat * m_R.col(1);
+    return T2Matrix() * m_R.col(1);
 }
 
 Eigen::Vector2cd MatrixRTCoefficients_v3::getKz() const
@@ -171,10 +156,10 @@ Eigen::Matrix2cd MatrixRTCoefficients_v3::computeInverseP() const
 Eigen::Matrix2cd MatrixRTCoefficients_v3::computeDeltaMatrix(double thickness)
 {
     Eigen::Matrix2cd result;
-    const auto alpha = 0.5 * thickness * (m_lambda(1) + m_lambda(0));
+    const complex_t alpha = 0.5 * thickness * (m_lambda(1) + m_lambda(0));
 
-    const auto exp2 = Eigen::Matrix2cd(Eigen::DiagonalMatrix<complex_t, 2>(
-        {GetImExponential(thickness * m_lambda(1)), GetImExponential(thickness * m_lambda(0))}));
+    const Eigen::Matrix2cd exp2 = Eigen::DiagonalMatrix<complex_t, 2>(
+        {GetImExponential(thickness * m_lambda(1)), GetImExponential(thickness * m_lambda(0))});
 
     // Compute resulting phase matrix according to exp(i p_m d_m) = exp1 * Q * exp2 * Q.adjoint();
     if (std::abs(m_b.mag() - 1.) < eps) {
@@ -182,15 +167,12 @@ Eigen::Matrix2cd MatrixRTCoefficients_v3::computeDeltaMatrix(double thickness)
         const double factor1 = 2. * (1. + m_b.z());
         Q << (1. + m_b.z()), (I * m_b.y() - m_b.x()), (m_b.x() + I * m_b.y()), (m_b.z() + 1.);
 
-        result = Q * exp2 * Q.adjoint() / factor1;
-
-    } else if (m_b.mag() < eps) {
-        result = Eigen::Matrix2cd::Identity() * GetImExponential(alpha);
+        return Q * exp2 * Q.adjoint() / factor1;
 
-    } else
-        throw std::runtime_error("Broken magnetic field vector");
+    } else if (m_b.mag() < eps)
+        return Eigen::Matrix2cd::Identity() * GetImExponential(alpha);
 
-    return result;
+    throw std::runtime_error("Broken magnetic field vector");
 }
 
 namespace
-- 
GitLab