Interface proposal

go.sh

@@ -2,7 +2,7 @@
 echo Compile with gcc -O2
 rm -rf *.bin
 
-# g++ -O2 -std=c++11 standalone.cc nmie.cc -lm -o scattnlay.bin
+g++ -O2 -std=c++11 standalone.cc nmie.cc -lm -o scattnlay.bin
 
 # clang++ -g -O1 -fsanitize=address  -fno-optimize-sibling-calls -fno-omit-frame-pointer -std=c++11 standalone.cc nmie.cc -lm -o scattnlay.bin
 

nmie-wrapper.cc

@@ -29,7 +29,6 @@
 ///
 /// @brief  Wrapper class around nMie function for ease of use
 ///
-#include "ucomplex.h"
 #include "nmie-wrapper.h"
 #include "nmie.h"
 #include <cstdio>

nmie-wrapper.h

@@ -1,4 +1,4 @@
-#ifndef SRC_NMIE_NMIE_WRAPPER_H_
+#ifNdef SRC_NMIE_NMIE_WRAPPER_H_
 #define SRC_NMIE_NMIE_WRAPPER_H_
 ///
 /// @file   nmie-wrapper.h
@@ -52,26 +52,97 @@
 
 namespace nmie {
   class MultiLayerMie {
+    // Will throw for any error!
    public:
+    // Set parameters in applied units 
     void SetWavelength(double wavelength) {wavelength_ = wavelength;};
-    void SetSizeParameter(std::vector<double> size_parameter);
-    void SetIndex(std::vector< std::complex<double> > index);
-    void RunMie(double *Qext, double *Qsca, double *Qabs, double *Qbk);
-    std::vector< std::vector<double> >  GetSpectra(double from_WL, double to_WL,
-                                                   int samples);
-    /// Disabled functions or unused untill the end of C++ conversion
+    // It is possible to set only a multilayer target to run calculaitons.
+    // For many runs it can be convenient to separate target and coating layers.
+    // Per layer
+    void AddTargetLayer(double layer_width, std::complex<double> layer_index);
+    void AddCoatingLayer(double layer_width, std::complex<double> layer_index);
+    // For all layers
+    void SetTargetWidth(std::vector<double> width);
+    void SetTargetIndex(std::vector< std::complex<double> > index);
+    void SetCoatingWidth(std::vector<double> width);
+    void SetCoatingIndex(std::vector< std::complex<double> > index);
+    void SetFieldPoints(std::vector< std::vector<double> > coords);
+
+    //Set parameters in size parameter units
+    void SetWidthSP(std::vector<double> width);
+    void SetIndexSP(std::vector< std::complex<double> > index);
+    void SetFieldPointsSP(std::vector< std::vector<double> > coords);
+
+    // Set common parameters
+    void SetAnglesForPattern(double from_angle, double to_angle, int samples);
+    std::vector<double> GetAngles();
+    
+    void ClearTarget();
+    void ClearCoating();
+    void ClearLayers();
+
+    // Applied units requests
     double GetTotalRadius();
+    double GetTargetRadius();
+    double GetCoatingWidth();
+    std::vector<double>                  GetTargetLayersWidth();
+    std::vector< std::complex<double> >  GetTargetLayersIndex();
+    std::vector<double>                  GetCoatingLayersWidth();
+    std::vector< std::complex<double> >  GetCoatingLayersIndex();
+    std::vector< std::vector<double> >   GetFieldPoints();
+    std::vector<std::vector<std::complex<double> > >  GetFieldE();
+    std::vector<std::vector<std::complex<double> > >  GetFieldH();
+    std::vector< std::vector<double> >   GetSpectra(double from_WL, double to_WL,
+                                                   int samples);
+    double GetRCSext();
+    double GetRCSsca();
+    double GetRCSabs();
+    double GetRCSbk();
+    std::vector<double> GetPatternEk();
+    std::vector<double> GetPatternHk();
+    std::vector<double> GetPatternUnpolarized();
+    
+
+
+    // Size parameter units
+    std::vector<double>                  GetLayerWidthSP();
+    // Same as to get target and coating index
+    std::vector< std::complex<double> >  GetLayerIndex();  
+    std::vector< std::vector<double> >   GetFieldPointsSP();
+    // Do we need normalize field to size parameter?
+    /* std::vector<std::vector<std::complex<double> > >  GetFieldESP(); */
+    /* std::vector<std::vector<std::complex<double> > >  GetFieldHSP(); */
+    std::vector< std::vector<double> >   GetSpectraSP(double from_SP, double to_SP,
+                                                   int samples);
+    double GetQext();
+    double GetQsca();
+    double GetQabs();
+    double GetQbk();
+    double GetQpr();
+    double GetAsymmetryFactor();
+    double GetAlbedo();
+    std::vector<std::complex<double> > GetS1();
+    std::vector<std::complex<double> > GetS2();
+    std::vector<double> GetPatternEkSP();
+    std::vector<double> GetPatternHkSP();
+    std::vector<double> GetPatternUnpolarizedSP();
+
+    void RunMieCalculations();
+    void RunFieldCalculations();
   private:
-    /// Size parameter for each layer.
-    std::vector<double> size_parameter_;
-    /// Complex index value for each layer.
-    std::vector< std::complex<double> > index_;
     const double PI=3.14159265358979323846;
-    /// Disabled functions or unused untill the end of C++ conversion
     void GenerateSizeParameter();
+    void GenerateIndex();
     double wavelength_ = 1.0;
     double total_radius_ = 0.0;
-
+    /// Width and index for each layer of the structure
+    std::vector<double> target_width_, coating_width_;
+    std::vector< std::complex<double> > target_index_, coating_index_;
+    /// Size parameters for all layers
+    std::vector<double> size_parameter_;
+    /// Complex index values for each layers.
+    std::vector< std::complex<double> > index_;
   };  // end of class MultiLayerMie
+
 }  // end of namespace nmie
 #endif  // SRC_NMIE_NMIE_WRAPPER_H_

nmie.cc

@@ -758,8 +758,11 @@ int nMie(int L, int pl, std::vector<double> x, std::vector<std::complex<double>
     *Qext += (n + n + 1)*(an[i].real() + bn[i].real());
     // Equation (28)
     *Qsca += (n + n + 1)*(an[i].real()*an[i].real() + an[i].imag()*an[i].imag() + bn[i].real()*bn[i].real() + bn[i].imag()*bn[i].imag());
-    // Equation (29)
-    // We must check carefully this equation. If we remove the typecast to double then the result changes. Which is the correct one??? Ovidio (2014/12/10)
+    // Equation (29) TODO We must check carefully this equation. If we
+    // remove the typecast to double then the result changes. Which is
+    // the correct one??? Ovidio (2014/12/10) With cast ratio will
+    // give double, without cast (n + n + 1)/(n*(n + 1)) will be
+    // rounded to integer. Tig (2015/02/24)
     *Qpr += ((n*(n + 2)/(n + 1))*((an[i]*std::conj(an[n]) + bn[i]*std::conj(bn[n])).real()) + ((double)(n + n + 1)/(n*(n + 1)))*(an[i]*std::conj(bn[i])).real());
     // Equation (33)
     Qbktmp = Qbktmp + (double)(n + n + 1)*(1 - 2*(n % 2))*(an[i]- bn[i]);